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Statistics of a General Class of Avalanche 
Detectors With Applications to Optical 
Communication 


By S. D. PERSONICK 
(Manuscript received June 29, 1971) 


Previous results on the statistics of avalanche detectors are generalized to 
the case where electrons and holes suffer collision tonizations with unequal 
probability. It is assumed here that the ratio of collision ionization prob- 
abilities per unit length of weaker-to-stronger carrier 1s a constant k in- 
dependent of position in the high-field region. The moment-generating 
function of the random avalanche gain G ts obtained as a function of k 
and the average gain G, and is used to obtain Chernov bounds on error rates 
of digital optical receivers employing avalanche detectors. It 1s shown the 
required energy per pulse to achieve a given error rate decreases as k de- 
creases for fixed G. For each k > 0, there is an optimal mean gain Guy. 
resulting in minimum required energy per pulse. Atk = 0.1, Goo. & 100 
and the required energy is within 10 dB of that required with very high 
gains (a few thousand) atk = 0. 


I. INTRODUCTION 


In a previous paper’ results on the statistics of two particular ava- 
lanche detectors with applications to optical communication were 
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presented. It was required either that only one carrier suffer collision 
ionizations in the high-field region (unilateral gain) or that both carriers 
suffer collision ionizations with equal probability per unit length in the 
high-field region. The present work allows for more general unequal 
ionization probabilities per unit length with the requirement only that 
the ratio of the two quantities be constant throughout the high-field 
region. The moment-generating function of the random gain is obtained 
as a function of this ratio and the average gain. The results are con- 
sistent with unpublished conjectures of R. J. McIntyre.” The moment- 
generating functions are used to obtain Chernov bounds on the error 
rates of digital optical receivers employing avalanche detectors and 
using either coherent or incoherent light. Results on avalanche statistics 
are summarized in Section V. Numerical results on the Chernov bounds 
are given in Section VI (6.5) and Section VII. 


II. MODEL OF THE AVALANCHE DETECTOR 


The avalanche detector is a device in which thermally or optically ~ 
generated hole-electron pairs generate additional hole-electron pairs 
through collision ionizations. Within the device there is a “high-field 
region” where holes have probability 8(x) per unit length (which depends 
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Fig. 1—Avalanche detector. 
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upon the position x) of suffering a collision ionization as they travel to 
the left under the influence of the electron field (see Fig. 1). Electrons 
traveling to the right have a probability a(x) per unit length of collision 
ionization. Carriers can be created within the high-field region due to 
thermal effects or due to the presence of incident light. Carriers can 
also drift into the region if they are generated outside of the region. It is 
assumed that all collision ionizations are independent. This requires 
that the mean distance between ionizing collisions be large compared 
to the distance over which a carrier can randomize its momentum after 
a collision.* Hole-electron pairs created through collision ionization can 
in turn generate additional pairs by the same mechanism. This is the 
avalanche process. 


Ill. THE STATISTICS 


We seek the statistics of the random total number of hole-electron 
pairs which result ultimately through collision ionizations when an 
initial hole-electron pair is injected into the high field at some position z. 
Define p,(n, x) as the probability that n pairs ultimately result in- 
cluding the initially injected pair. The moment-generating function of 
the number of pairs M,(s) is therefore* 


M,(s, 2) = Diplo, ae”. (1) 


We shall derive M,(s, x). Before proceeding we must review some well 
known results which will be needed in that derivation. 

If {x,;} are random variables which are independent, then the moment- 
generating function of the sum of the {xz,} is the product of the individual 
moment-generating functions.* The semi-invariant moment-generating 
function SIMGF of a random variable X having probability density 
px(x) is the natural logarithm of the moment-generating function of X 


¥x(s) = In [Mx(s)] = In b pia" | : (2) 


The SIMGF of a sum of n independent random variables is the sum of 
the individual SIMGF’s. 

We can now proceed to derive M,(s, x). Divide the high-field region 
into K intervals of width dX = W/K. See Fig. 2. Label these intervals 
1, 2, 3, --- j, --: K. If a hole-electron pair is injected into interval j, 
define, as above, the probability density of the total number of pairs 
ultimately resulting in the avalanche process (including the initial pair) 
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as p,(n, x) where x is taken as the center of interval 7. The hole of the 
initial pair moves to the left and the electron to the right toward 
x = Oorz = W respectively. As they pass through their respective 
intervals, new pairs may be created in each through collision ionization. 
We shall assume that the interval width dX is sufficiently narrow so 
that the initial pair carriers create either one or no new pairs in each 
interval. If the initial hole or electron generates a new pair in some 
interval k, then that new pair will ultimately generate NV, pairs including 
itself through the avalanche process. Thus with each interval we can 
associate a number of pairs N,. This number equals zero if the appro- 
priate initial pair carrier suffers no collision ionizations in interval k. 
This number equals one or more if the appropriate initial pair carrier 
suffers a collision ionization in interval k. The total number of pairs 
ultimately generated through the avalanche process including the 
initial pair is one plus the sum of the {N,}. Since collision ionizations 
are all independent, all the N, are independent. Thus we have the 
SIMGF of the total number of pairs given by 


W/dx 
¥,(s, x) 25 + p> vw, (s) (3) 
where s is the SIMGF of the deterministic initial pair and wy,(s) is 
the SIMGF of N,. 
The SIMGF of NV, is obtained as follows. The probability that. 
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Fig. 2—Avalanche process. 
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N, = Oisl1 — y dX wherey = a(z) if interval k is to the right of interval 
j where the initial pair enters, y = B(z) if k is to the left of j. (z is eval- 
uated at the center of interval k.) The probability that NV, = Z > Ois 
the probability of a collision ionization in interval k times the probability 
that a new pair created at interval k ultimately results in Z pairs 
including itself. That is 


Pr (N, = Z) = 7 dX p,(Z, x) for Z>0 (4) 


with x evaluated at the center of interval k. Thus 
Wy,(s) = In bs Pr (N,. = ze | 
Z=0 


=In la —ydX)e° + >) y dX p,(Z, ae? | 
Z=1 


=Infl-—ydX + y7dX M,(s, x)] 


=Infl —ydX +ydXe""”] (5) 
where 
y=a(z) if k>j 
y= Bz) if k <j 


x = center of interval k. 


Using (8) and (5) and taking the limit as dX gets infinitely small, one 
obtains 


vs, 2) = s+ | Ble’ fe? — 1 da’ 
: Ww , 
+ fo aleve"? — 1) de’. ©) 
Equation (6) is the critical equation for determining y,(s, x), and thus 
M,(s, x) = exp [y,(s, x)]. Using Leibnitz’s rule for differentiation of 


integrals one obtains 


= dal, 2) = (82) — alae" — 1). (7) 


The solution of (7) is 


¥.(6, 2) = In 7 @) 
1 — C exp ( [B(2’) — a(z’)] dz’) 
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where C = (e¥9'? — 1)/e¥*"*" which can be checked by substitution. 
Substituting (8) into (6) one obtains the particular result 


WY, (s, 0) 
1 


=s+ i a(x’) anc na (J ae es 2 — 1| dz’. (9) 


If one makes the asswmption® that at each point in the high-field 
region 
B(x) = k-a(zx) (10) 


where k is a constant, one can solve (9) to obtain 


1 gre 
g\?) 0 95 sec | | 
° ° y k—1 M,{s, 0) > eM, (s, 0) - 1] 


= 8§ + ar In [M,(s, 0) = eT Ms, 0) — 1]] (11) 
where 
Ww 
M,(s,0) =e"? and 6= / a(x) dx. 
0 


One can write (11) in another way by making a substitution.* Define 
6’(s) implicitly by 


ee = M,(s, Oe’, (12) 
that is, 
ev) — Mi(s,0) = e PF, (13) 
Using (12) in (11) one obtains the implicit equation 
me ee A GLa Pi en a (14) 


which determines 6’(s) and thus M,,(s, 0) through (13). 

Equation (14) is still not explicit. A numerical technique for solution 
is discussed in the next section. One can use (11) [or (14)] and (8) to 
obtain y,(s, z) or M,(s, x) for any zx. Recall that x is the point of entry 
of the initial pair. In the applications we shall be concerned with + = 0 
or x = W. That is, pairs are generated in a drift region outside the high- 
field region with carriers drifting into the high-field region. 

* Equation (14) will follow from (11) and (12) by tedious algebra. Further, (14) 


will not be used in the following results except to compare with McIntyre’s work 
in Appendix A. 
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The equation (14) is consistent with some unpublished conjectures of 
MelIntyre given in Appendix A. 


IV. NUMERICAL SOLUTIONS 


Equations (11) through(14) can be solved numerically. One technique 
is to differentiate (11) to obtain the result . 


= Ms, 0) = Mls, ober — § (Ms, ee ar (15) 


where V/,(0, 0) = 1. 
Equation (15) can be integrated with a computer to obtain /,,(s, 0) 
explicitly. 


V. SUMMARY OF ANALYTIC RESULTS ON AVALANCHE STATISTICS 


From Sections I through IV and a previous paper by this author,’ 
we obtained the following: 


Assumptions: 


Holes travel toward x = 0, electrons toward x = W. 
Hole ionization probability per unit length = B(x). 
Electron ionization probability per unit length = a(z). 
B(x) = k-a(x), k a constant for all z. 

High-field region width = W. 


Definitions: 


5= [ a(x) dx, 


p,(n, x) = probability that if an initial pair enters the high-field region 
at point x, n pairs will ultimately result through the avalanche process 
including the initial pair, 


G(x) = mean avalanche gain = >> np,(n, 2), 
1 


M,(s, x) = moment-generating function of p,(n, 2) 


> pal, ze 


Results: 
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LER E _ ye E +i) [ aw) at || 


forall k. (16) 
For k = 0 (Unilateral Gain) 


2a. M,(s, 0) = [1 — e*f1 — e']]7’ where G(0) = e’. (17) 
ob. 560 Oe 1 (ay where @ = GO). (18) 


Fork ¥0,k #1 


ca nat k —w 1 = 1 —s 5\k-1 re 
3. as M,(\s, 0) _— M,(s, 0( k Jf k (M,(s, O)e € | ’ 


where 
M,(0,0) = 1 (19) 


co (EN fom 


For k = 1 (Equal Ionization) 


where 


4. a M;(s, 0) = M(s, afi = (2 a, 0) | 
M,(0, 0) = 1 
where 
G = G0). (20) 


VI. APPLICATIONS TO RECEIVERS USING AVALANCHE DETECTORS 


6.1 General Comments 


We shall next apply the results of Section V to obtain bounds on the 
error rates of digital receivers. The receivers to be discussed here are 
the single- and twin-channel systems described below. We shall upper- 
bound the power required at the receiver to obtain a desired error rate 
' using the Chernov bounds. 


6.2 The Receivers 


The twin-channel receiver is shown in Fig. 3. Depending on the state 
of a binary information source, one of two channels has optical output 
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Fig. 3—Twin-channel system. 


power p for duration T while the other has output power p-E XT for 
duration 7. EXT signifies extinction. EXT would ideally be zero, but 
is left finite and less than unity for practical reasons. The optical power 
falling on.an avalanche photo diode causes the emission. of photo- 
electrons which are multiplied along with the detector ‘dark current”’ 
through the avalanche gain mechanism. The detector outputs (multi- 
plied counts) are integrated in devices having thermal noises referred to 
their respective inputs. The integrator outputs are subtracted and the 
difference X is compared to a threshold of zero to decide what the 
information state was. 

The single-channel system is essentially half the twin-channel system 
as shown in Fig. 4. The single integrator output X is compared to a 
threshold y to decide upon the information state. 


6.3 The Chernov Bounds 


The Chernov bound is a useful tool for bounding the probability that 
a random variable X will lie above or below a given threshold y. It is 
given as follows” 


Pr BS > 7) < el X(s)- sw x(a) Dade os 


provided s > 0, 
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UxGISavet ) , 
Pr (X < 7) Se wpe 


provided s < 0, (21) 


where 


¥x(s) = SIMGF of X = In i px(x)e” as | 


and 
ve) = S ox). 


6.4 Chernov Bounds for the Two Systems 


6.4.1 Preliminaries 


We wish to determine the required power p to achieve a desired 
decision error probability for each of the receivers discussed in Section 
6.2 with various types of detectors and various values of other system 
parameters such as dark current and thermal noise in the integrators. 
To apply the Chernov bounds we will need the SIMGF’s of the variables 
X at the outputs of the receivers. (See Figs. 3 and 4.) 

An important result needed here is the following: 


Lemma: If C ts tnteger-valued and greater than or equal to zero; and if 
U = eS gi where the g; are independent, identically distributed random 
variables (that is, each ‘‘count’’ produced by the C process independently 
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Fig. 4—Single-channel system. 
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generates g; contributions to the U process through some gain mechanism); 
then the SIMGF of U is’ 


Wu(s) = vcly,(s)) (22) 
where w,(s) ts the SIMGF of the tdentically distributed random variables g; . 


One can model an avalanche photodiode as a photon counter followed 
by a random avalanche multiplier. Each photon counter output ‘‘count”’ 
produces a random number of “counts” at the multiplier output. Thus, 
since we know from Section V the SIMGF of a random multiplier which 
corresponds to y,(s) in (22), and since we seek the SIMGF of the photo- 
diode output which corresponds to Wy(s) in (22), it follows that we need 
the SIMGF of the number of counts emitted by a photon counter with 
light incident upon it, which corresponds to W¢(s) in (22). 


6.4.2 Photon Counter Statistics 


6.4.2.1 Coherent Light. If the light incident upon the photon counter 
is coherent, then it is well known that the SIMGF of the total counts 
emitted in an interval 7’ is given by® (see Appendix B) 


Wc(s) = [LAMSIG + LAMD] [e’ — 1] (23) 
where 
LAM SIG = Total incident light energy -7/hQ 
n = Detector quantum efficiency 
hQ = Energy of photon at optical frequency used 


LAMD Mean number of dark current counts before avalanche 
gain in an interval T. 


6.4.2.2 Incoherent Light. If the incident light is incoherent with H 
independent spatial-temporal “degrees of freedom,” then the SIMGF of 
the total counts emitted by the photon counter in interval T is (see 
Appendix B) 


Wo = LAMD[e’ — 1] + In IE — aes (e — »}"| (24) 


where LAMSIG is the average total incident energy times 7/hQ. 
Clearly (24) is the same as (23) as H approaches infinity, which is a 
well known result. 


6.4.3 Final Calculations 


6.4.3.1 T’win-Channel System. Since, for the twin-channel receiver, 
X consists of the difference of two random integrator outputs, we need 
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the following well known result.* If X = X, — X,, and if X, and X, 
are independent, then the SIMGF of X is 


¥x(s) = px,(s) + ¥x,(—8). (25) 


Each integrator output contains the sum of the counts emitted by its 
detector and the integral of its thermal noise. The SIMGF of the random 
variable N obtained when Gaussian thermal noise spectral height N, is 
integrated over an interval T is well known to be 


2 
by(s) = pl: (26) 
Using (22) through (26) we obtain for the SIMGF of the twin-channel 


receiver output X when the information is in state ‘fone’, and the 
optical source is coherent, 





¥x(s) = NT + [ern BE LAMD |1a1.(9 —q 


a |eE-EXT 9 os LamD |[M(—9 —1) (27) 


hQ 
where 
LAMD = mean number of dark current counts before avalanche 
gain in an interval T 
N, = spectral height of the thermal noises referred to the 


integrator inputs. 


M,(s) is obtained from (16) through (20) depending upon the particular 
gain mechanism. If the optical source is incoherent, we have 


¥x(s) = 9° N oT + LAMD[M,(s) + M,(—s) — 2] 


+ In IE — 211M, @) — nf "| 
+in|{ [1-27 2X? yy — 0) "]- es 


We seek the probability that when the information is in state ‘‘one’’, 
X is less than zero, and we therefore decide that the information state 
was ‘“‘zero.”’ That is, we seek the error probability. One can use (27) 
or (28) and the Chernov bound of (21) to determine the required value 
of LAMSIG = p-T-n/(hQ) to achieve a desired error probability. 
Since the twin-channel receiver is symmetric, the error probability 
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when the information is in state “zero” is the same as when it is in 
state ‘“‘one.”’ 


6.4.3.2 Single-Channel System. Tor the single-channel receiver, we 
need the SIMGF of X under both information states. Call X, the random 
variable X when the information is in the state ‘‘one.”’ Call X, the 
random variable X when the information is in state zero. Using results 
of Section 6.4, one obtains for coherent light 


vee) = NE | namo +221 Nor - 1 


ve) = Me 4 | ramp +P EXP aay — 1. es) 


For incoherent light 


Net S Nol oe 7 AM DIM G1] 


+ In [1-2 — aan 
Net 


-)"] 
ves) = AS + LAMDIM,() — 1] 


+ In IE — 2 Ee (M6) nf "|. (30) 


One can then use the results of (29) and (30) along with the Chernov 
bounds of (21) to simultaneously find values of LAM SIG = p-T-7/(hQ) 
and the threshold y (see Fig. 4) to ensure some desired error probability 
(which for convenience here will be the same for either information state). 


vx.) = 








6.5 Numerical Results 


The Chernov bounds described above were evaluated numerically. 
The results are displayed on the attached figures described below. The 
range of parameter values is realistic and practical, to the best of this 
author’s knowledge. The curves presented are those deemed most 
interesting by the author. Other calculations can of course be made. 
Parameters used are defined as follows:* 


LAMSIG = Required mean number of detected photons per 
pulse in the “on” channel of the twin-channel 


* SIG, EXT, G, K, H, and LAMD are input parameters to the program which 
calculates LAM SIG for a desired error rate. 
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receiver or in the ‘‘one”’ state of the single-channel 
receiver. 

LAMSIG.EXT = Mean number of detected counts per pulse in the 
“off?’? channel of the twin-channel receiver or in 
the ‘‘zero”’ state of the single-channel receiver. 

SIG = Normalized thermal noise standard deviation 


= {4R0T/[Re’]}* = {440C/e"}* 


where e = electron charge, 6 = Boltzmann’s 
constant-absolute temperature, R = equivalent 
noise resistance at integrator input, 7’ = pulse 
duration, C = T/R = integrator equivalent input 
capacitance. For the results to follow, a reasonable 
value of SIG was chosen to be 6000. 
= Mean avalanche gain. 
= Temporal-spatial diversity for incoherent carrier 
case. 
k = Ratio of ionization probability per unit length of 
weaker and stronger ionizing carriers.* 
LAMD = Dark current counts per interval T' before 
avalanche gain. 


Md 
| 


Fig. 6 


LAMSIG vs G is plotted for the twin-channel case with k as param- 
eter. SIG was set at 6000, the error rate is 10°, LAMD was set at 
5 counts and EXT = 0.01. H = 10,000 which is equivalent to assuming 
a coherent carrier. 


Fig. 6 


The value at optimal gain of LAMSIG vs k is plotted. Points are 
tagged with the optimal G. The receiver is a twin-channel system with 
SIG = 6000, EXT = 0.01, LAMD = 5. H is 10,000 which is equivalent 
to assuming a coherent carrier. The error rate is 107°. 


Fig. 7 
LAMSIG vs G is plotted for two values of error rate 10°° and 10°° for 
* For these calculations it was assumed that the detector is designed so that the 
stronger ionizing carriers generated optically or associated with dark current enter 


the high-field region from a drift region outside the high-field region. This corresponds 
to initial pairs entering the gain mechanism of x = 0 orz = W as discussed in Section 
III. 


LAMSIG 


STATISTICS OF AVALANCHE DETECTORS 3089 


6000 


SIG = 6000 
5000 LAMD = 5 
EXT = 0.01 
H = 104 (COHERENT) 
TWIN CHANNEL 
ERROR RATE = 107° 





4000 
[e) 
a 
= 3000 
< 
a 
2000 
1000 
re) 
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Fig. 6—LAMSIG at optimal gain versus k. 
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Fig. 7—LAMSIG versus gain. 


a twin-channel system with SJG = 6000, EXT = 0.01, H = 10,000 
LAMD = 5,k = 0.1. 
Fig. 8 


Single- and twin-channel systems are compared. LAMSIG vs G is 
plotted for SIG = 6000, EXT = 0.01, LAMD = 5, H = 10,000, error 
rate = 10°’, k = 0. Note that from an average power viewpoint the 
single-channel system is 3 dB better than shown if the binary informa- 
tion source is random, since LAMSIG is the energy in “‘one”’ state. 


Fig. 9 

Same as Fig. 8 except k = 1. Note the scale change. 
Fig. 10 

LAMSIG vsG for H = 100 and H = 10,000 for twin-channel system. 
SIG = 6000, EXT = 0.01, LAMD = 5, k = 0, error rate = 10°”. 
6.6 Further Comments 


When systems were investigated for sensitivity to the choice 
LAMD = 5, EXT = 0.01, it was found that insignificant changes in 
LAMSIG vs G occurred when various combinations of LAMD = 5 
or 50, EXT = 0.01 or 0.001 were tried. 
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Fig. 8—LAMSIG versus gain. 
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Fig. 9—LAMSIG versus gain. 
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VII. CONCLUSIONS ON APPLICATIONS 


If one assumes that the Chernov bounds are sufficiently tight so that 
actual energy required per bit to achieve specified error rates can be 
compared for various system parameters by comparing the bounds,* 
then one can conclude the following. 


(i) Define k as the ratio of the collision ionization probabilities per 
unit length of the weaker-ionizing to the stronger-ionizing carrier 
(carriers are of course holes and electrons). Assume that the detector 
is designed ‘‘well’’ such that optically and thermally generated carriers 
enter the high-field region from a drift region outside. From the bounds, 
one obtains the result that the required energy per pulse to achieve a 
desired error rate decreases as k decreases for fixed average avalanche 
gain. A value k = 01s best; but a value & = 0.1 will allow one to operate 
with energy within 10 dB of that required at very high gains with a 
k = 0 device. For each value of k except zero there is an optimal gain 
resulting in minimum required energy per pulse. The optimal gain is 
larger for smaller k. At k = 0.1, the optimal gain is about 100. At 
k = 0, the optimal gain is infinite, but a gain of a few thousand allows 


* For simple cases where both the bounds and actual energy requirements can 
be obtained (for instance for the k = 0 case) the two results differ by a few dB or 
less. 
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close to optimal required energy per pulse. One can conclude that a 
silicon device with k = 0.1 and a gain of about 100 would be a good 
choice for an optical detector. This is true since a detector with k less 
than 0.1 and yet having gain significantly higher than 100 is not avail- 
able at this time. 

(11) The required energy per pulse for systems using incoherent 
optical sources differs from that for systems using coherent sources by 
less than a few dB provided the product of the source bandwidth and 
the pulse duration exceeds 100. This is true even if there is no spatial 
incoherence of the light at the detector. 

(217) For reasonable parameter values, and assuming a random 
information stream, the single-channel receiver requires about 1.5 dB 
less energy per pulse to achieve a desired error rate than the twin-channel 
receiver. 

(iv) The required energy per pulse is insensitive to reasonable values 
of dark current and extinction ratios. 

(v) For a particular system, a change in the desired error rate from 
10°° to 10°° results in a change in the required energy per pulse of 
1 to 3 dB, depending upon the avalanche gain. This shows that the 
required energy per pulse is fairly insensitive to the error rate. On the 
other hand, this means that poor error rates will result if insufficient 
loss margin is provided. That is, a small lowering of the received energy 
can greatly increase the error rate. 


APPENDIX A 


In an unpublished work, McIntyre conjectures (from special case 
calculations) that the probability density of the random gain, defined 
here as p,(n, 0) is given by 





r(; i i a. ere _ a a 


! u 2—n) 
nite +2 a 





p(n, 0) = 





where k and 6 are the same as in (10) through (13). 

If one makes the assumption that the conjectured p,(n, 0) has sum 
over 7 normalized to unity for each value of k and for each 6, then one 
obtains the result of (14) by using the definition of the moment-generat- 
ing function and the normalization property. 
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APPENDIX B 


If light of known intensity falls upon a photon counter during an 
interval 7’, then the probability density of the total number of counts 
emitted is well known’ to be Poisson distributed as follows 

go erTh) 
pe(n) = [A + LAMD' —_—. (31) 

Where A* is the total energy incident in the interval T times 7/AQ, 
LAMD is the mean number of dark current counts per second times 
the interval 7, and 7/AQ is the detector quantum efficiency divided by 
the energy in a photon. . 

The moment-generating function of the distribution of (31) is given by 


M.(s) = exp [[A + LAMD]le’ — 1]]. (32) 


If the incident light is a stochastic process, then the moment-generat- 
ing function of the output count distribution is obtained by averaging 
(82) over the probability density of the stochastic total energy incident 
in the interval T 


Mes) = [exp (A + LAMDIle" — 1)p(a) aa. (33) 


An incoherent light field is normally taken to mean that the complex 
envelope of the classical field is a complex Gaussian random process. 
That is, such a field incident on the photon counter plane can be written 
as 


E(p, t) = v2 refe(o, te" } 
pe counter plane (34) 
te (0, 7) 


where e¢(p, ¢) is a complex Gaussian random process. 
If one expands ¢(p, ¢) in its Karhunen—Loeve eigenfunctions,’ one 
obtains 


e(p, i= >» CPx (Pp, t) 
pe counter plane (35) 


te (0, 7’) 


* In the text, A is called LAMSIG. 
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where 


T 
[ fF aslo, D6%6, 0 ap at = a, 
counter plane “0 


and the coefficients e, are independent complex Gaussian random vari- 
ables satisfying 


(exes) = Vebk.,; (36) 
(e,e;) = 0. 
The energy incident upon the photon counter is 
12 n= | elo, De, Dd pdt = ¥ lel? (37) 
: p, te*(p, p €;,| « 


If one assumes an equal distribution of average energy in roughly H 
‘“‘modes,”’ 


= < < 
Yk Y) l1sksdH (38) 


= 0, k>H 


then it follows that from (83) and the complex Gaussian statistics of 
the e, that 


M-.(s) = exp [LAMD@€ — y-[t — 5a ver — yf" (39) 


Assumption (38) implies that the energy of the incoherent light is 
roughly equally distributed in H degrees of freedom. 
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Design of Digital Filter Banks for Speech 
Analysis 


By R. W. SCHAFER and L. R. RABINER 
(Manuscript received June 22, 1971) 


A bank of bandpass filters 1s often used in performing short-time spectrum 
analysis of speech signals. This paper 1s concerned with the analysis and 
design of digital filter banks composed of equally spaced bandpass filters. 
It is shown that significant tmprovement in the composite filter bank 
response can be achieved by proper choice of the relative phases of the band- 
pass filters. The results are extended to more general filter bank configura- 
tions. 


I, INTRODUCTION 


Many speech processing systems are based on the concept of short- 
time spectrum analysis.’’” Spectrum analyzers for such systems often 
consist of a set of bandpass filters whose combined passbands cover a 
desired frequency range. Although continuous-time filters have tradi- 
tionally been used in filter banks for speech analysis, hardware realiza- 
tions of digital filters are now available,’ and the advantages which 
digital filters offer should be exploited in filter bank design. These 
advantages include: flexibility of design of the individual bandpass 
filters, precision of realization, stability of digital hardware, and the 
efficiency of realization of the filter bank afforded by the possibility 
of multiplexing the digital hardware. Thus it is important to consider 
design techniques for filter banks composed of digital bandpass filters. 

To focus on the basic concepts in filter bank design, it is useful to 
define an zdeal filter bank spectrum analyzer. Figure 1 depicts such a 
filter bank composed of digital filters whose impulse responses are 
denoted by h,(nT), k = 0, 1, --- , M, where 1/T is the sampling 
frequency of the input signal.* Such a filter bank constitutes an ideal 
spectrum analyzer if the input 2(n7) (with possibly further band limit- 
ing) can be synthesized exactly (within some fixed delay) by a linear 


* The filter ho(n7') is a lowpass filter which is included for completeness although 
this band is usually not analyzed in practical speech analysis systems. 
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Fig. 1—Bank of digital bandpass filters. 


combination of the bandpass filter outputs 2,(n7). An example of 
such a system would be one in which the filters are ideal rectangular 
bandpass filters with the same constant gain and linear phase in their 
passbands and zero gain outside. If the filter bandwidths are such 
that the frequency range 


—~T 


i 


is completely covered without overlap, then the input can be synthesized 
exactly by adding together the outputs of the bandpass filters. 

The essential characteristics of the ideal spectrum analyzer are that 
the frequency response of the combined outputs must exhibit a flat 
magnitude response and a linear phase, and therefore the combined 
impulse response must be a delayed digital impulse (unit sample). 
Causal digital filters (filters whose impulse responses are zero for 
nm < 0) cannot have the desired ideal gain characteristics and may 
not have linear phase.* Therefore a filter bank composed of such filters 
cannot achieve the ideal characteristics of flat magnitude response and 
linear phase. In this paper we will describe an approach to the design 
of filter banks that approximate the ideal spectrum analyzer. First 
we present a detailed analysis of a filter bank configuration in which 
equally spaced, equal-bandwidth digital filters are used. This analysis 
suggests a technique for optimizing the filter bank characteristics and 


lA 


Tv 
e=T7 


* Finite duration impulse response digital filters can have precisely linear phase. 


DIGITAL FILTER BANKS 3099 


also suggests how the results can be used where nonuniform bandwidths 
are desired. The results are illustrated with examples. 
II. ANALYSIS OF UNIFORM FILTER BANKS 


Assume that the bandpass filters in Fig. 1 have impulse responses 
of the form 


h,nT) = 2 | D, | h(nT) cos @nT + &,) k=1,2,---,M 
ho(nT) = | Dy | h(n) 


where w, = w, + (kK — 1) Aw, and h(nT) is the impulse response of a 
prototype lowpass filter. (When the k = 0 lowpass filter is used, w, = Aw 
and w) = 0.) The system functions for this set of bandpass filters are 


H(z) = | D, | e’®*H(ee***”) + | D, | eo 7?* A (ze'***) k=1,2,---,M 
H(z) = | D. | H@). (2) 


The frequency response of these filters is obtained after substituting 
z = e'*” in (2) as 


Hylel**) = | Dy | el HEP?) + | Dy [eH EY, 
k=1,2,---,M, (3) 


(1) 


H,(e'*") = | Do | H¢’*"). 


If the frequency response of prototype lowpass filter, H(e’°”), drops 
off sharply, then it can be seen from (8) that 


le") | ~~ | D; | LHe 2) | 0 < a < «/T 
w|D,||HE°**?*") |—«/T So SO. 


In this case the filter bank consists of a set of (J7 + 1) equally spaced 
bandpass filters with identical magnitude responses around their 
respective center frequencies. We have chosen this method of designing 
bandpass filters from lowpass prototypes for analytical convenience 
and because of the importance of spectrum analysis systems of this 
form.’"” The results to be discussed apply for other bandpass trans- 
formations in so far as they yield a set of uniformly spaced bandpass 
filters with identical frequency characteristics. 

Our objective is to choose the prototype lowpass filter and the 
parameters | D, |, w,, Aw, and &, so that the filter bank will closely 
approximate the characteristics of the ideal spectrum analyzer. To do 
this we must consider the response of the composite filter bank. First, 
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however, it is useful to interpret the individual bandpass filter outputs 
in terms of spectrum analysis considerations. 
The individual filter outputs are of the form 


2,(nT) = be 2e(rT) | Dy | h(n — rP) cos [wx(nT — rT) + 4), (4) 


Rn) = > BET) De hae) 


r=—-o 


which can be expressed as 


z,(nT) = 2 Re {D,X@, ,né°"7} 


(5) 
x(nT) = Re {D.X(0, n)} 
where D, is the complex constant defined by 
D, = | D, | e** k=1,2,-+-,M (6) 
Do = | Do |, 
and 
X@,n) = DS arT)AMT — rT)?” (7) 


The quantity X(«,, ”) is the discrete-time version of the short-time 
Fourier transform* of z(n7). Thus, (5) serves to relate the bandpass 
filter outputs x,(n7T) to the short-time Fourier transform. 

The frequency response and impulse response of the composite filter 
bank are obtained from 


M 


y(nT) = >» x,(nT). (8) 


k=0 


After substituting (5) into (8) and noting that if x(n T) is real, X(—w, , n) 
is the complex conjugate of X(w, , n), we obtain 


M 
y(nT) = Dd) DX ,n)ei**"" (9) 
k=—-M 
where w_, = —w,, and D_,is the complex conjugate of D, . Substituting 


(7) into (9) and interchanging the order of summations results in 


n 


y(nT) = »> x(rT) Kaa te = Denor |. (10) 


Defining 
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AM 
d(nT) = a>) D,e***, (11) 


we observe from (10) that the combined impulse response of the filter 
bank can be expressed as 


A(nT) = h(nT) d(nT). (12) 


Equations (11) and (12) are the basic results of the analysis. (Note 
that they could have been obtained directly by summing the impulse 
responses h,(n7T), with the sacrifice of the interpretation of the filter 
bank outputs in terms of the short-time spectrum.) 

Equation (12) shows that the filter bank impulse response is the 
product of the prototype lowpass filter impulse response h(nT), and 
the sequence d(nT’) defined by (11). The choice of lowpass filter depends 
on both the desired frequency resolution and the requirement of obtain- 
ing flat magnitude and linear phase response in the composite filter 
bank. The sequence d(nT’) is independent of the prototype lowpass filter 
and is a function of the frequency spacing, the relative gains and phases, 
and the number of bandpass filters. Thus, for a given choice of prototype 
lowpass filter, the parameters of d(nT) can be adjusted to achieve the 
best approximation to the ideal spectrum analyzer. To see how this 
occurs, we shall first examine in detail the characteristics of the sequence 
d(nT). 

As will be shown in the remainder of this section, a particularly useful 
choice of the complex coefficients D,, in (11) is 


D, =e?" k= +1,42,---,4M. 
Do = 1, 


(13) 


where 7 is an integer. That 1s, 
®, = @,NoT’, and | D, | = 1. 


(The condition Do = 1 implies that the band around w = 0 is included 
in the filter bank; D) = 0 implies that it is not.) It can be shown that 
for this choice of D, , (11) becomes 


2sin [M Aw(nT + noT’)/2] 
sin [Aw(nT' + n,T’)/2| 





d(nT) = Do + cos [w,(n7' + noT')| (14) 
where w, = w, + (M — 1) Aw/2, and D, is 1 or 0 depending on whether 
or not the lowpass channel is included. 

The properties of the sequence d(nT) determine the character of 
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the impulse response of the filter bank. Some of these properties are 
summarized below: 


(t) The parameter > shifts the sequence d(n7’) by no samples with 
respect to h(nT). 


(ii) The sequence d(nT) is even about the sample n = —m;i.e., 
d(nT — nT) = d(—nT — nT). 
(iit) The maximum value of d(nT) occurs at nm = —n,. This value 


is d(—nmT) = 2M+ D.. 

(iv) If 2r/(AwT) and w,/Aw are both integers, then the sequence 
d(nT) is periodic with period 27/Aw. Otherwise, d(nT) will be 
an almost periodic sequence which will peak up at time intervals 
of 27/ Aw. 


Insight into the properties of d(n7T’) can be gained by considering a 
simple example. Assume that w, = Aw and 27/(Aw7’) = N where N is 
an odd integer. That is, the entire frequency range —7/T Sw S 2/T, 
is divided into N equal bands. If 14 = (N — 1)/2, the entire frequency 
range is covered. Under these conditions (14) can be written 


i | 2) Aw(nT + net) | 
sin [Aw(nT + nT)/2] 


sin | 2H) enh net | 


sin [Aw(nT + noT’)/2] 


(If N is even, the k = 0 filter is not used, i.e., Do = 0, and w, is Aw/2.) 
It is clear from (15b) that for these conditions d(n7) is a periodic 
sequence with period NT = 27/Aw. In fact, d(nT) may be thought of as 
samples of a continuous-time periodic Dirichlet kernel as shown in 
Fig. 2a. If M = (N — 1)/2 and D, = 1, d(nT) is a periodic discrete- 
time impulse train with impulses occurring at multiples of N 7. This 
is because the sample points on the periodic Dirichlet kernel occur at 
the maxima and the zero crossings, as indicated by the small circles in 
Fig. 2a. 

The conditions for d(nT) to be periodic are that both w,/Aw and 
27/(AwT’) be equal to integers. To see this, we must examine (14) in 
detail. If 2r/(Aw7') is an integer, and M is an odd integer, the sequence 
2 sin [M Aw(nT + noT)/2]/sin [Aw(nT + noT)/2] is periodic with 
period NT = 22/Aw. If w,/ Aw is an integer and M is odd, the sequence 
cos [(w, + (M — 1) Aw/2)(nT + np T)] is also periodic with a period 


Il 


d(nT) = —-1 if D=0 (15a) 


if Do =1. (15b) 
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Fig. 2—(a) Periodic continuous-time Dirichlet kernel, (b) continuous-time en- 
velope and sequence d(n7') when either w/ Aw or 27/(Aw7') are not integers. 


that is an integer multiple of 2r/Aw. Thus the product of these two 
sequences is periodic with period 27/Aw. The identical result holds 
for M an even integer and 27/Aw and w,/Aw integers although the 
interaction between the component sequences is slightly different. 

If either 27/(AwT7') or w,/Aw are not integers, d(nT) will not be 
periodic, but will still peak up at time intervals of 27/Aw. Such a case 
is depicted in Fig. 2b where the samples d(nJ) are marked by the 
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small circles and the dotted curve shows the factor 
sin [M Aw(nT + noT")/2)/sin [Aw(nT + no T)/2] 


when M is odd. As shown in Fig. 2b, d(nT) will always have even sym- 
metry about sample n = —n. 


III. DESIGN OF UNIFORM FILTER BANKS USING BESSEL FILTERS 


In the preceding section we presented a detailed analysis of a filter 
bank composed of equally spaced equal-bandwidth filters. In this 
section we will show how the results of that analysis can be employed 
in filter bank design. 

The objective of flat amplitude response and linear phase is most 
easily achieved with bandpass filters having these same properties. 
For this reason, Bessel (maximally flat delay) filters are often used in 
filter banks.’ In the examples shown in this paper, we have used digital 
filters obtained from Bessel prototype designs using impulse invariance.° 
It should be noted that the digital filters obtained this way do not have 
the maximally flat delay property. J. P. Thiran’ has shown that the de- 
nominator of the system function of maximally flat delay digital filters 
is a Gauss hypergeometric function. It is reasonable to expect however, 
that for the narrow-band filters of interest here, the differences should 
be negligible. 

As an example a digital filter derived from a sixth-order Bessel 
lowpass filter with asymptotic cutoff frequency of 60 Hz is shown 
in Fig. 3. The impulse response is shown in Fig. 3a, and the amplitude 
and phase responses are shown in Fig. 3b and Fig. 3c. The filter shown in 
Fig. 3 was used in a filter bank* with the following choice of parameters: 
D, = 0, T = 10% sec, Aw = 27(100), w, = 27(100), nm = 0, and 
M = 30. The resulting filter bank characteristics are shown in Fig. 4. 
The filter bank impulse response, h(nJ), is shown in Fig. 4a along 
with the prototype lowpass impulse response h(n7'). For the above 
choice of parameters, d(n7’) is obtained from (15a) as 


sin [0.6177] _ 


AO) tole) (16) 
which is periodic with period 100 samples (10 msec), with peaks occurring 
at nT = 0, +10, +20, --- msec. From Fig. 4a, it can be seen that 


in the product h(nT)-d(nT), the peak of d(nT) at nT = O will be 
attenuated since h(nT) is small around nT = O. On the other hand, 


* Note that the resulting bandpass filters are twelfth order. 
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_Fig. 3—Sixth-order Bessel filter characteristics. (a) Impulse response, (b) mag- 
nitude response, (c) phase response. 
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Fig. 4—Characteristics of 30-channel filter bank. (a) Impulse response (dotted 
curve is the impulse response of the prototype lowpass filter in Fig. 3), (b) composite 
magnitude response, (c) composite phase response after subtracting 10-msec delay. 
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the peak of d(nT) at nT’ = 10 msec occurs at approximately the peak 
of h(nT), and at nT = 20 msec, h(nT) is large enough to produce a 
significant echo in the impulse response of the filter bank. As is shown 
in Fig. 4b and 4c, this corresponds to a 3.9-dB ripple in the amplitude 
response and a 25.5-degree peak-to-peak ripple in the phase response 
(after removing a linear phase component corresponding to a 10-msec 
or 100-sample delay). To decrease this amplitude and phase ripple, we 
should attempt to eliminate the echo in the impulse response. Further- 
more, the phase ripple will be eliminated if the impulse response h(nT) 
has even symmetry about some delay time np7'. One approach is to 
broaden the filter bandwidths, or equivalently reduce the spacing Aw, 
so that h(nT) is contracted relative to the spacing of pulses in d(nT). 
This is generally not an acceptable solution since h(nT) and Aw are 
usually fixed by some frequency resolution criterion. However, if we 
refer to the properties of d(nT) which were previously summarized, 
we note that a negative value of no will shift d(nT) to the right relative 
to h(nT) so that d(nT) will have even symmetry about time npT = 
—oT'. If np can be chosen so than h(nT) = h(nT)-d(nT) has approxi- 
mately even symmetry and consists of only one significant pulse, then 
the amplitude and phase ripple will be small. The manner in which this 
is achieved is shown in Fig. 5 where it is assumed for simplicity that 
d(nT) is a train of digital impulses as would be the case for M = (N — 
1)/2. Figure 5a depicts the case where n) = 0. Figure 5b shows the 
situation where n) was chosen to shift the impulse which was at nT’ = 0 
in Fig. 5a to the right and into the vicinity of the peak of h(nT). If it 
is assumed that only three impulses have nonzero amplitudes (a, , 
Q2 , a3) such that 4 | a |-| as | < | a, + a |-| a, |, then it can be shown 
(see Appendix) that the peak-to-peak amplitude ripple of the filter 
bank is 


R, = 20 logy, eect |], (17) 


Similarly, if | a, + a3 | < ||, the peak-to-peak phase ripple about 
a linear phase corresponding to a delay of —7 7’ is given by 


=1 a, — A 
ges la ea + rl ai 
The conditions for (17) and (18) to hold are satisfied when a, and a; 
are small relative to a, , which is the normal situation. It can be seen 
from (18) and (17) that the phase ripple will be zero if a, = a3, and 
the amplitude ripple will be small if (a, + a3)/a, is small. 

Although these results were derived for the idealized case when d(nT) 
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Fig. 5—Illustration of how to adjust the parameter no. (a) Composite impulse 
response for no = 0, (b) no chosen to minimize magnitude and phase ripple (dotted 
lines indicate movement of individual pulses in d(n7')). 


is an impulse train, we have found that amplitude and phase ripple 
can be determined quite accurately using (17) and (18) in more general 
situations. With the foregoing principles in mind we have written an 
interactive computer program for filter bank design. Using this program 
we can design a filter bank with low amplitude and phase ripple by 
the following process: . 


(t) Choose w, , Aw, and M to cover the desired analysis band and 
choose h(nT) to provide desired frequency resolution. This 
results in an h(nT) that has a duration of approximately 471/Aw 
as shown in Fig. 5. | 

(ii) Evaluate h(nT) and determine mn) such that a, & a3; as in 
Fig. 5b. . 

(iz) If the resulting filter bank is not satisfactory, steps 7 and 7 
are repeated. 


In cases where w,/Aw is not an integer, it is important to choose no 
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so that the point of even symmetry in d(nT’) is shifted into the vicinity 
of the peak of h(nT’). Otherwise, it may be impossible to achieve a very 
good approximation to linear phase. An example of the improvement 
gained by proper choice of nm) is shown in Fig. 6. In this example all 
the parameters were the same as in the example of Fig. 4 except a value 
of ny = —129 was chosen by the above process. In this case the in-band 
amplitude ripple is 0.8 dB and the phase ripple is 0.6 degree, as compared 
to 3.8 dB and 25.5 degrees when n) = 0. 

R. M. Golden’ has shown that inverting the sign of alternating chan- 
nels often significantly improves the characteristics of a filter bank. This 
technique has a simple interpretation in terms of our results. It can be 
shown that inverting the sign of alternating channels is equivalent to 
delaying the sequence d(nT) by no = —7/(AwT’) samples. This amount 
of delay may be nearly correct if the duration of h(nT) is approxi- 
mately 32/Aw; however for the situation shown in Fig. 4a, such a delay 
would produce a worse filter bank than no delay at all (n> = 0). Also, 
to achieve linear phase when w,/Aw is not an integer, the point of even 
symmetry in d(nTJ) should be delayed to the vicinity of the peak of 
h(nT). This does not occur when the signs of alternate channels are 
inverted. 


IV. DESIGN OF NONUNIFORM BANDWIDTH FILTER BANKS 


In speech applications it is common to take advantage of the fre- 
quency resolution characteristics of the ear*’® by using increasing 
bandwidth filters at higher frequencies. The previously discussed tech- 
niques can be applied to this situation if the filter bank consists of 
several sub-banks, each with different resolution. Each sub-bank can 
be designed as discussed above, with care being taken to ensure that 
the entire frequency band of interest is covered by the combination 
of the sub-banks. It may be necessary to equalize the delay between 
sub-banks by providing additional delay for all but one of the sub- 
banks.* This is depicted in Fig. 7 for three sub-banks with increasing- 
bandwidth sixth-order Bessel filters. Figure 7a shows the lowpass 
prototype impulse response and shifted d(n7T) sequence’ for the first 
sub-bank. The lowpass asymptotic cutoff used was 78 Hz, the spacing 
of filters was Aw, = 27(125), the first filter was centered at w,, = 
27(250), and a value of nm), = —100 (10-msee delay) was required to 


* Golden’ has shown that the delays can be approximately equalized by increasing 
the order of the lowpass prototype in direct proportion to the increase in bandwidth. 

* The sequence d(n7’) is shown as an impulse train for convenience in plotting. 
The actual sequences would look like those in Fig. 4 and Fig. 6. 
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Fig. 6—Characteristics of 30-channel filter bank. (a) Impulse response for no = 
—129 (dotted curve is the impulse response in Fig. 3), (b) composite magnitude 
response, (c) composite phase response after subtracting 12.9-msec delay. 
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minimize the amplitude and phase ripple. Figure 7b shows the second 
sub-bank in which the basic parameters were: lowpass asymptotic 
cutoff 136 Hz, Aw. = 27(218), wi. = 27(1296.5), and m2. = —57 
(5.7-msec delay). To line up the central peaks, an additional delay 
of n. = 43 samples (4.3 msec) was required. Figure 7c shows the third 
sub-bank where the lowpass cutoff was 192 Hz, Aw; = 27(807), wi3 = 
27(2213), and m3; = —40 (4-msec delay). A value of n; = 60 samples 
(6.0 msec) is required to line up the central peak with those in Fig. 7a 
and 7b. The response of the combination of these three sub-banks is 
shown in Fig. 8. Figure 8a shows the impulse response, Fig. 8b shows 
the amplitude response, and Fig. 8c shows the phase after a linear 
phase corresponding to 10-msec delay has been subtracted. It can be 







nT (msec) 


nT (msec) 





® 10 20 nT(msec) 


Fig. 7—Illustration of the design of nonuniform filter banks: (a) impulse response 
for narrow bandwidth filters, (b) impulse response for intermediate bandwidth filters, 
(c) impulse response for wide bandwidth filters. 
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_ Fig. 8—Composite filter bank characteristics with three different sub-banks: (a) 
impulse response, (b) magnitude response, (c) phase response after subtracting 
10-msec delay. 
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seen in Fig. 8b and 8c that the ripple in the sub-banks is quite low 
as would be expected from Fig. 7. At the boundary between sub-banks, 
however, the ripple increases significantly due to the fact that the last 
filter in the lower sub-bank drops off more rapidly than the first filter in 
the next sub-bank. This excessive variation at the boundary between 
sub-banks can be eliminated to some extent by using increasingly 
higher-order filters in the sub-banks. Alternatively, nonuniform resolu- 
tion can be obtained by using equal-bandwidth filters and adding 
together groups of two or more of their outputs to achieve the desired 
bandwidth. Such an approach would require increased computation 
but would produce filter bank characteristics comparable to those in 
Tig. 6. 


Vv. CONCLUSION 


We have discussed the analysis and design of digital filter banks 
and have shown how the incorporation of a linearly increasing phase 
shift in each bandpass filter can significantly improve the overall filter 
bank characteristics. We also showed how the techniques can be used in 
nonuniform bandwidth filter banks. 

The examples which we gave were based on Bessel lowpass prototypes 
which have impulse responses of desirable shape but rather poor am- 
plitude response. Recent results in the design of finite duration impulse 
response filters® offer attractive possibilities for filter bank design. 
Such filters can have precisely linear phase and can be designed using 
iterative techniques with constraints on both the impulse response 
shape and the amplitude response. The use of such filters, together with 
the basic principles discussed in this paper, should yield filter banks 
with excellent properties. 


APPENDIX 
Derivation of Magnitude and Phase Ripple Formulas 
Assume an impulse response sequence 
h(n) = a n=0 

= a, n= Ny 

= 3 n = 2n, 

= 0 elsewhere. (19) 
The system function of this system is 


H(e?**) = a, + ae ??"?? + aye i0??”, (20) 
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The squared magnitude response is 
| H(e'*7) |? = fas + (a, + a) cos (wn, T)? + (ax — a3)” sin” (wn,T), 
(21) 


and the phase response is 


arg [H(ei*")] = tan” eee | ea 


where a linear phase component —wn,7 has been removed. Clearly, 
both (21) and (22) are periodic functions of w with period 27/n,T. 
To determine the amplitude and phase ripple, we must locate the 
maxima and minima of (21) and (22). 

If we differentiate (21) with respect to w, we find that the maxima 
and minima occur for values of w satisfying 


sin (wn,T') = 0 (23a) 
(a, + a3) a Os) (23b) 


COs (wn, 1") = —az ea 
1 


The second equation is satisfied by a real value of w if and only if 
4] a |-[as | > | or + ae |-| a |. (24) 


In a good filter bank design, a, and a; will be positive and much smaller 
than a, , and (24) will not be satisfied. Evaluating the second derivative 
shows that in this case the maxima and minima of | H(e’*”) | will 
alternate and occur at values of w satisfying (23a); i.e., 0 = 0, -72/n,T, 


+27r/n,T, ---. In this case the amplitude ripple in dB is given by 
R, = 20 logy, [faete tal], (25) 
~ Ll ae — ay — as | 


If (22) is differentiated with respect to w, we find that the maxima 
and minima occur at values of w satisfying 


cos wn,l’ = — (ata), (26) 
a2 

Equation (26) is satisfied by real values of w if |a: + a3| < | a |. 

In this case the maxima and minima again alternate, and the peak- 

to-peak phase ripple is 


ms a Qa, — a3 oe 
Rem Bean [ee | a 
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If |a, + a3| > | a.|, the phase curve will be discontinuous with a 
jump of 27 radians occurring at w = +a/n,T, +30/n,T, --- . 
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The Preference of Slope Overload to 
Granularity in the Delta Modulation of 
Speech 


By N. S. JAYANT and A. E. ROSENBERG 
(Manuscript received June 18, 1971) 


A preference study was made to assess the relative annoyance values 
of slope-overload distortion and granular noise in delta-modulated speech. 
A recently described adaptive delta modulator was simulated at frequencies 
of 20 and 40 kHz, and controlled amounts of the two types of degradation 
were introduced into samples of a 2-second utterance. Rankings were 
obtained for these samples on the basis of preference judgments of nine 
listeners, each of whom assessed the samples, pairwise, in a tournament- 
type strategy. Results indicate that the speech sample exhibiting the mini- 
mum degradation on an objective, overall-notse-power basis 1s not subjectively 
the most preferred sample. Furthermore, the subjectively optimum delta 
modulator exhibits greater overload and lesser granularity than the ob- 
jectively optimum device. 


I. INTRODUCTION 


The principle of delta modulation’ has been widely described in the 
literature. Briefly, delta modulation is a digital encoding strategy 
which uses a simple feedback mechanism to produce a “staircase” 
approximation to an input signal. A block diagram of the simplest form 
of delta modulation appears in Fig. 1. The input sequence {X,} is 
usually band-limited and suitably oversampled. The ‘‘staircase”’ se- 
quence Y, is generated according to the equations 


C, = sgn (X, — Y,-;) (1) 
Ye eg ee Me SS By SO (2) 


The step-size A, is assumed to be a constant in conventional (linear) 
delta modulation. ‘‘Adaptive” delta modulation, on the other hand, 
allows for modifications of A, in accordance with the changing slope 
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Tr= Xp Yp-1 
\ pgn Trp= Cp 
/ 
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Fig. 1—Schematic diagram of a linear delta modulator. 


characteristics of the input signal. Such adaptation results in better 
encoding, and several types of adaptive delta modulation have been 
described in the literature.”'*’* 

Figure 2 illustrates the mechanism of an adaptive delta modulator 
and demonstrates how suitable increases and decreases of step size 
facilitate better encoding during steep and flat regions of the input 
signal waveform. Such adaptations can be effected by observations on a 
“recent”? segment of the binary sequence {C,}; this is illustrated by 
equation (5) in the sequel. 

Figure 2 also brings out the distinction between two types of encoding 
error in delta modulation, viz., “granular noise” and ‘‘slope-overload’”’ 


"HUNTING" OR "GRANULAR" 
NOISE 


"SLOPE-OVERLOAD" 
DISTORTION 
aN 





Xp- Xp_y=Sp 






Ae Ye Ete: Xt-17 Yt-1 7 E t-vo 
sgn [Xp- Yr-1]=Cp= sgn [mr] 
Cp= TRANSMITTED "CHANNEL” SYMBOL 


Fig. 2—Illustration of adaptive delta modulation. 
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distortion. A given error sample 
E,= X, —Y, (3) 


can be defined to fall into the granular or slope-overload category, 
depending on whether the corresponding step m, crosses the input 
waveform or not. Thus, in Fig. 2, there is a ‘granular’ error Eg at the 
sampling instant ¢, and an ‘overload’ error E(,-,)9 at the sampling 
instant (¢ — 1). As a matter of definition, we will note that H,, = 
Bau-ne = 9. 

The signal output {Z,} of the delta modulator is actually obtained by 
filtering the staircase sequence {Y,} to the input signal band. Let 
{X"} be the result of passing {X,} through the same lowpass filter. A 
perceptually relevant measure of signal degradation is accordingly 
defined by the encoding error 


t= X= 7. (4) 


As with the quantity #, in (3), one can distinguish samples of granularity 
and slope overload, e,g and e,o, in the error sequence {e,}. Referring 
to Fig. 2 once more it can be seen that a physical distinction between 
the two types of error is suggested. Granularity can be described as a 
“signal-uncorrelated’’ random noise-type of phenomenon. It is char- 
acterized by alternation of signs and tends to be independent of signal 
amplitude. Slope overload, on the other hand, can be described as a 
“signal-correlated”’ distortion, since its sign and magnitude are related 
to the slope of the signal. This physical difference between slope overload 
and granularity suggests a corresponding perceptual distinction and 
raises the question of the relative annoyance values of the two forms 
of signal degradation in delta modulation. The present paper describes 
a study of the above question as referred to the delta modulation of a 
speech signal. 

Earlier work in this subject is in the form of a perceptual experiment’ 
in which H. Levitt, et al., characterized the perceptibility of slope- 
overload distortion as such. As mentioned earlier, our paper will seek 
to answer the complementary question of the relative perceptibilities 
of slope overload and granularity when they occur simultaneously in 
delta-modulated speech, as they usually do. 

The approach we used was to vary the relative amount of slope 
overload and granularity introduced into samples of a test utterance, 
and to evaluate these samples on the basis of both objective and per- 
ceptual criteria; and then to interpret these evaluations with specific 
reference to the overload-granularity dichotomy. 


3120 THE BELL SYSTEM TECHNICAL JOURNAL, DECEMBER 1971 


Section IJ summarizes the salient features of a computer-simulated 
adaptive delta modulator that was utilized in the present study. This 
adaptive encoder has been recently described and shown to provide 
toll-quality speech reproduction at bit rates of practical importance.* 

Section III defines the objective measures of speech quality used in 
our study, while Section IV defines a subjective measure of preference 
and describes an underlying perceptual experiment. 


II. DESCRIPTION OF THE DELTA MODULATOR 


Figure 3 is a schematic block diagram of the adaptive delta modulator 
utilized in the present study. This encoder is defined by the basic 
equations (1) and (2), and by the adaptation rule 


Ava Pr Ay TC. = C,& 
1 ; -P 


= p Ari if C, ¥ Cri 


IV 
— 


(5) 


Notice that a conventional (linear) delta modulator corresponds to 
the special case of P = 1. In our study the value of P was a variable 
parameter; different (delta-modulated) speech samples corresponded to 
different suitably spaced values of P, and thereby to different mixtures 
of slope-overload and granularity. 

The original speech sample X was a 2-second male utterance of ‘‘Have 
you seen Bill?” that had been band-limited to 3.38 kHz. The delta 
modulation was performed at sampling rates of 20 and 40 kHz. The 
latter frequency provides speech reproduction that approaches telephone 
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Fig. 3—Schematic diagram of an adaptive delta modulator. 
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quality.* The lower sampling rate was included to provide a better 
demonstration of the annoyance properties of delta-modulated speech. 


III. OBJECTIVE MEASURES OF SPEECH QUALITY 


We recall the encoding error e, , (4), and define the following measures 
of delta-modulator performance. (Summations are over the entire length 
of the speech utterance, and a nonzero granularity error at r = t implies 
zero overload error, and vice versa.) 


(7) The overload-noise energy in Z: 


No = oD CG | (6) 
(21) The granular-noise energy in Z: 
Nes > es (7) 


(iit) The signal-to-noise ratio: 





pae.¢ 
SNR = 3 2 8 
>» C0 + > Org ( 


(iv) The signal-to-granular (overload)-noise ratio: 
xX; 


> eat 


IV. A SUBJECTIVE MEASURE OF SPEECH QUALITY 


SNRe co) = 





(9) 


The perceptual evaluations. of this paper are based on the pooled* 
judgments of nine listeners each of whom assessed speech stimulit 
in six runs of a perceptual experiment. Each of these 54 experiments 
was a double-elimination tournament? (with a different, random, 
starting line-up). Matches in each tournament were between contending 
stimuli, playing two at a time. The result of each match was in the 
form of a binary preference judgment by the listener, while the result 
of a tournament was a set of scores awarded to each of the contesting 
speech stimuli on the basis of its record in the tournament. The actual 
scoring rule’ was one which, together with the double-elimination 

_* Intralistener variations were found to be less significant than the intrastimulus 
differences. 

t The number of contending speech stimuli was also nine, at each sampling rate. 

+ The tournament ended when every losing contestant had lost twice. 

§ In the course of each tournament a contestant accumulated a score as follows. 
No. score was earned for a match that was lost; while, after every match that was 


won, the contestant’s score was the sum of the accumulated scores, before the match, 
of the contestant and of the loser, plus one. 
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strategy, provided a useful alternative—as concluded from a separate 
simulation—to the more comprehensive testing procedure where every 
contending stimulus would be pitted against every other. 

It was recognized, however, that both the scoring rule and the double- 
elimination strategy were empirical procedures. This was more so 
because they were applied to what was apparently a probabilistic 
environment: the binary preference-response of a listener to a given pair 
of contending stimuli can well be random, especially when the stimuli 
are not obviously different. It was, therefore, decided not to emphasize 
the actual scores obtained in the perceptual test. They were only used, 
instead, to extract a crude ranking information that would be less 
sensitive to the testing and scoring procedures. 

Consequently, the following subjective preference value Q was 
assigned to each of M contesting speech stimuli: 


M—-R 
li 7 R=1,2,---,M (10) 


where RF is the rank assigned to a stimulus on the basis of its accumulated 
score in the 54 runs of the perceptual test. 


V. SUMMARY OF RESULTS 


Figure 4 displays normalized values of the objective measures of 
quality SNR, SNRz , and SNRo , as well as the subjective preference 
function Q, as functions of the adaptation parameter P. The following 
observations emerge: 


(c) The speech sample representing the minimum overall-noise-energy 
is not subjectively the most preferred sample. In fact, at both 20 and 
40 Hz, the objective and subjective optima can be characterized by 


Poon = 12 (11) 

oon = 1b (12) 

(22) The approximate coincidence of the SNR and SNR» curves 
indicates, by virtue of equations (6) through (9), that 

No>WNe (13) 


for all considered values of P. 

(417) The relative disposition of the SNR», SNR¢, and Q curves— 
and of their maxima—demonstrates that, in spite of the preponderance 
(13) of overload in the overall-noise-energy, the granularity in a speech 
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NORMALIZED PERFORMANCE 





Fig. 4—Evaluations of delta-modulator performance. 


stimulus has a strong influence on its subjective preference value. 

Note that there is a double peak in the subjective preference curve, 
Q, for the 40-kHz case. This curve unambiguously ranks each of the 
experimental stimuli according to equation (10). However, the actual 
scores underlying this ranking show only a small difference between 
the stimulus with the secondary peak and the one immediately preceding 
it. What is probably indicated is a general broadening of the peak of 
the preference function for values of P between 1.2 and 1.5. 

Table I lists, for the optimal characterizations (11) and (12), values 
of Noa (given as fractions of signal energy), SNR, and Q. Notice that 
the subjectively optimum delta modulator displays lesser granularity 
(Ng) and greater overload (N 0) than the objectively optimum modulator. 
It is again obvious that in perception, overload and granularity are 
not weighed in proportion to the respective noise energies No and Ng ; 


3124 THE BELL SYSTEM TECHNICAL JOURNAL, DECEMBER 1971 


TABLE I—CHARACTERISTIC OF OPTIMAL ADAPTIVE DELTA MopuLATION 
(No and Ng are entered as fractions of signal energy) 








Sampling 
Frequency P No Ne SNR Q. 
SUBJ 
P = 1.2! 0.0216 0.0003 43 1 
OPT 
20 kHz 
OBJ 
P =1.5 | 0.0158 0.0004 58 0.81 
OPT 
SUBJ 
P = 1.2] 0.0022 0.00003 450 1 
OPT 
40 kHz 
OBJ 
P = 1.5} 0.0016 0.00004 640 0.91 
OPT 


in fact, the perceptual preference of a speech sample seems to be deter- 
mined very strongly by the extent of granularity in it, although the 
latter represents a very small fraction of the total noise energy. — 

. Finally, Table I indicates that distinctions between objective and 
subjective assessments of speech quality appear to be less significant 
at the higher sampling rate of 40 kHz; thus, for example, the objectively 
best delta modulator has a greater value of subjective perference Q 
at 40 kHz than at 20 kHz. 


VI. CONCLUSION 


We have shown that in delta modulation, a speech sample exhibiting 
the minimum degradation on an objective, overall-noise-energy basis 
is not equivalent, in general, to the perceptually most preferred sample. 
We have also indicated that this distinction may be less significant in 
higher quality delta modulation than in a low-bit-rate encoder. 

The subjectively optimum delta-encoder displays a greater overload 
No and lesser granularity Ng than the objectively best encoder. This 
feature, together with the fact that No > Ng in either case, suggests 
the strong influence of granular noise on the perceptual assessment 
of a speech sample; equivalently, a lesser ‘‘annoyance value”’ is to be 
associated with slope-overload distortion.* A possible explanation of 


* Companding in PCM exploits a similar but not identical subjective phenomenon, 
viz., the greater tolerance to encoding errors in regions of high input amplitude. 
(Notice however, that in delta modulation, slope overload is not confined to high- 
amplitude regions, nor is granularity associated only with low input amplitude.) 
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this observation would be the fact that granularity is explicitly per- 
ceivable by a listener as an ‘‘additive background noise,”’ while slope- 
overload distortion exists only in relation to an original signal which is 
not known to the listener. 

Finally, our observation that slope overload is “‘less annoying”’ than 
granularity is to be invoked with caution. Broadly speaking, we believe 
that our conclusion would apply very well to speech that achieves or 
approaches telephone quality. In extremely low-quality delta modula- 
tion (such as may be used in special applications), on the other hand, 
the intelligibility of speech will be a critical criterion; and in such a 
situation, depending on other factors like ambient noise at a transmitter, 
slope overload may very well become a more important perceptual 
attribute. 
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Some Considerations of Error Bounds 
in Digital Systems 


By V. K. PRABHU 
(Manuscript received June 22, 1971) 


Simple upper and lower bounds on the distribution function of the 
sum of two random variables are presented in terms of the marginal distribu- 
tion functions of the variables. These bounds are then used to obtain upper 
and lower bounds to the error probability of a coherent digital system in the 
presence of intersymbol interference and additive gaussian noise. The 
bounds are expressed in terms of the error probability obtained with a 
finite pulse train, and the bounds to the marginal distribution function of 
the residual pulse train. Since the difference between the upper and lower 
bounds can be shown to be a monotonically decreasing function of the 
number of pulses in the finite pulse train, the bounds can be used to compute 
the error probability of the system with arbitrarily small error. Also when 
the system performance is evaluated by simulation techniques, the methods 
presented in our paper can be utilized to estimate the error caused by using 
a finite pulse train approximation. 


I. INTRODUCTION 


In digital transmission systems the transfer characteristics of the 
transmitting and receiving filters are far from ideal, and the real trans- 
mission channel usually exhibits some form of time dispersion.’’” When 
an ideal digital signal is passed through such filters or is transmitted 
through such a channel, the successive pulses overlap; this form of 
distortion is usually known as intersymbol interference. Intersymbol 
interference may also result from the choice of nonoptimum sampling 
instants, imperfect demodulating-carrier phase, improper pulse design, 
etc. In addition the signal may be corrupted by thermal noise, co- 
channel and adjacent channel interference, and other forms of noise 
that may be present in the channel or in the system used to transmit 
the information. 

In digital transmission systems, one of the main performance char- 
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acteristics is the probability of error; this probability of error can often 
be expressed as a finite weighted sum of one or more distribution 
functions. 

Various authors have tried to evaluate this probability of error by a 
variety of methods,” ‘* but this highly complex probability distribution 
can seldom be exactly computed. 

Simulation techniques that may be used to solve this and other similar 
problems are never exact since one is constrained to use only a finite 
number of pulses and no bounds to the truncation error have been 
derived.* 

Another method is an analysis by means of a worst-case or ‘‘eye 
pattern” analysis. Since the probability of occurrence of a worst sequence 
may be very small, this analysis usually leads to very pessimistic results 
and suboptimum system design. 

Recently, some authors have derived’ several different upper 
bounds on the probability of error when the system is subject to both 
intersymbol interference and additive gaussian noise. Some of these 
bounds make use of the Chernoff inequality in their derivation, and 
hence are often more useful than the worst-case bound.*® However, 
since these bounds, in certain cases, can be shown to be loose,’’ and 
since no useful lower bounds have been derived, they are not as useful 
in system design as the evaluation of the exact error rate of the system. 

The third method consists in using the finite pulse train approxima- 
tion and calculating the error probability either by the direct enumera- 
tion of all possible sequences” or by the series expansion method.’””™* 
The series expansion method, which involves the computation of the 
moments of the intersymbol interference, is a convenient method but 
is still inexact as no truncation error bounds due to the residual pulse 
train have been derived. Note that in this method the number of terms 
in the finite pulse train is gradually increased until the change in 
probability of error is less than a given number e.” 

In this paper we first present simple upper and lower bounds to the 
distribution function of the sum of two random variables zy and Zz 
in terms of their marginal distribution functions. If the spread or dis- 
persion’” of the random variable zz is smaller than the spread of the 


ies 


* In simulation techniques the number N of pulses are usually chosen so that the 
computed probability of error stops changing by less than e when the number N is 
increased by 1. Noting that the series )):~ 1/n diverges, and that the difference 
between two successive partial sums of this series can be made less than any given 
number e, one concludes that this technique of choosing N is mathematically un- 
sound. 
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random variable zy , one can show that these two bounds are fairly 
close to each other and that one can evaluate the distribution function 
of the sum of the variables in terms of the distribution function of zy 
and the bounds on the distribution function of z,» . 

We then use these bounds to obtain upper and lower bounds on the 
error probability of a binary coherent digital system in the presence 
of intersymbol interference and additive gaussian noise. Since the 
difference between the upper and lower bounds can be shown to be a 
monotone decreasing function of the number N of pulses in the finite 
pulse train, the bounds can be used to compute the error probability 
of the system with arbitrarily small error. 

Also when the system performance is evaluated by simulation tech- 
niques, the methods presented in our paper can be utilized to estimate 
the error caused by using a finite pulse train approximation. 

If the symbols are equally likely, we also show that another set of 
upper and lower bounds can be derived for the probability of error of a 
system subject to intersymbol interference and additive gaussian noise. 

The usefulness of the bounds is illustrated by two examples. 


II. DISTRIBUTION FUNCTION AND ITS EVALUATION 


Let us assume that a random variable z is the sum of two random 
variables zy and Zp , 


2=2y + 2p, (1) 
and that we are interested in the distribution function of z2 
F(a) = Pr [2 S a] = Pr [ey + Ze S al. (2) 


In this section we shall also assume that zy and 2p are statistically 
independent random variables. 

The probability of error of a large number of digital systems subject 
to various forms of noise can often be expressed as a weighted sum of 
F,(a)’s. If z is the sum of an infinite number of random variables, and 
if zy represents its partial sum of the first NV terms, we sometimes can 
evaluate F,,(a), but F,(a) can seldom be computed exactly. In such a 
case it is often advantageous to obtain upper and lower bounds to 
F(a) in terms of F,,(a) and some known parameters associated with 
the random variable zz , the sum of the remaining terms in z. If the 
difference between the two bounds is a strictly monotone-decreasing 
function of NV, we can then calculate F,(a) with arbitrarily small error. 

Without loss of generality we shall assume that the mean of 2, is 
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zero. From (2) we can write (see Fig. 1) 
PQ) = ff fensealos 9) de dy ) 
if the joint probability density function f,,.,,(z, y) exists; and 
F=f [ aPa@ Fu) = [ Pala- dF),  & 
or 


F(a) = (F. x(a = Jee : (5) 


Let us now select an interval (— Al, Au) from the range of the random 
variable zz . From (4) we can write 


Fia)=1I,+1,+ 1, (6) 
where 
—Al 
he f Fule- dF), G 
Au 
I, = - F,,(a = y) dF, ,(y), (8) 
-Al 
and. 
I= | F.s(a— vy) dF a). @) 
Au 
ZR 


> 


Fig. 1—Distribution function of z = zy + 2p. 


ERROR BOUNDS IN DIGITAL SYSTEMS 3131 
One can show (see Fig. 2) that 


—-Al 
oss] aF..Q) = F..(-A), (10) 


051, 5 F.,(a— Au) [ dF.) = Paya — Aed{1 — F,,(40)} 
= F,y(a ae Al){l o> F,,(Au)}, (11) 


I, & Fay(a— Au) [ dP.g(y) = Fag(a — 0) {F.,(Au) — Paola), 
a: 12 


and 


Ie S Fag(a + AD [ aP.n(Q) = Fen(a + AD{R.q(Au) ~ Feq(—AD)}. 


(13) 
Combining (6) with (10)-(13), we have 
F,,(a — Au)[F,,(Au) — F,,(— Al] 
S F,@) Ss F,,(—Al) + F.,(a + AD[I — F,,(—AD] 
S F,,(—Al) + F,,(a + Al). (14) 


In general it is not easy to compute F,,,(y). However we may be able 
to bound F,,(y) so that 


0 < F,,(—Al) = Pr [zr S$ —Al] € L,,(—Al) 


IIA 
- 


(15) 


051 -— F,,(Au) = Pr [ze > Au] S U,,(Au) (16) 


IIA 
ae 


and 
1 = F,,(Au) — F,,(—AlD = Pr[—Al < ze S Au] 
21-1L,,(—Al) — U,,(Au)20. = (17) 
If these bounds can be found, (14)—(17) can be made to yield 
F,y(@ — Au)[1 — L,,(—Al) — U.,(Au)] S F.(a) 
= PGP AD haf AD. ~ (18) 


These are the basic bounds that we shall use in the rest of this paper. 
If the mass of the distribution of zz is very much concentrated 
around y = 0, our technique of computing F(a) from (18) relies on 
the assumption that we can find two numbers Au and Al such that 
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Fz, (2) 
1 





(b) 


Fig. 2a—Distribution function Fz, (2x). 
Fig. 2b—Distribution Tinetion F.,(y). The interval (—Al, Au) is contained in 
the range of zz, and for all practical purposes the mass of zp is contained in (— Al, 
Au). 


Au «|a|, Al<|a|, L,,(—Al) « F,,(a), U.,(Au) K F,,(a), and 
F,,y(a — Au) + F,,(a + Al). 

The difference D(Au, Al) between the upper and lower bounds can 
be written as 


D(Au, Al) = Pr [a — Au < zy S a+ Alj 
is F,y(a = Au)[Lz 2(— Al) ae U,,(Au)] +" Ta(—A)- (19) 


If Aw and Al can be so chosen that they are strictly monotone-decreasing 
functions of N, Au — 0, Al > 0, as N > o, and if the bounds on the 
distribution of zz are such that, for sufficiently large N, L,,(— Al) 
and U,,(Au) can be made smaller than any given number e¢, , we can 
estimate F’,(a) from (18) with arbitrarily small error.* 

For any given N even though Aw and Al can be chosen by optimizing 
the bounds in (18), this optimization leads to very complex equations. 
Hence we think that an algorithm should be developed to choose Au 
and Al for any given zy and zz. The development of this algorithm 
will be illustrated by an example in Section IV. 


* We assume that Pr[a —0 < ay Sa+0] =0. 
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2.1 Lower Bound Evaluation with Convex F(a) 

We shall now derive a simpler lower (upper) bound to F,(a) if F,,,(a) 
is a convex (concave) function and if zp is an even random variable, 
or 


FA(-t) =1— F,@. (20) 


From (20) one can show that the mean m of z is zero, and that its 
probability density f,(¢), if it exists, satisfies the equation - 


f.(—t) = f(t). (21) 


If 2 is an even random variable, we shall set Aw = Al in (18). 

Let us now assume that zp 1s an even random variable and that F’,,,(a) 
is convex over the range (a — Au, a + Au) where (— Au, Aw) is the 
range of zz . Since Zz is an even random variable 


F(@) = (Fin(@ — Yer)ben = (Fen@ + Yee) den (22) 
or 
F(a) = (3[Fvla — yea) + Pen(@ + Yse))ee - (23) 
Since F’,, (a) is convex over the range (a — Au, a+ Au),’° 
PF v(@ — Ys2) + Fiyla + y.2)] & Fen (Q)- (24) 
From (23) and (24) we have 
F(a) 2 F,,y(a). (25) 


Since this bound does not contain Au and Al, it is simpler to calculate 
than that given in (18). It is also tighter than the lower bound in (18). 
In this case we then have 


F,y(a) S$ F.(a) S$ F.,@ + Al) + L,,(— Al). (26) 


If F.,(a) is concave over the domain (a — Al, a + Au) and if zp 
is an even random variable, we can similarly show that 


F,,(@) 2 F,(@) 2 F,y(@ — Au)[l — L,,(— Al) — U,,(Au)]. (27) 


2.2 Evaluation of Another Upper Bound to F ,(a) 


Often we find that 2 contains a gaussian random variable n and can 
be written as 


z=n+ wy + ze = ent ee, zn=n+un, (28) 


where n, wy, and Zz are statistically independent random variables. 
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We have already assumed that the mean of zp is zero. Without loss 
of generality we shall now assume that the mean of n is zero, and its 
variance is o°. 

From (28) one can show that 


F,,(a) = 1 erte (tte), (29) 
where | 
aieays a | ee (Pak (30) 


Hence we have’”’* 


F,,(a) = 4 erfe (=) + Fr exp [—a@’/20"] 


| ps (—1)'H,_,(—a/oV3)(1/oV3)'m,/k!, (81) 


where H,(x) is the kth order Hermite polynomial and py, is the kth 
moment of wy , 


m= | dP, (2). (32) 
If the range (—Q, , 2,) of wy is finite and if Q denotes the maximum 


absolute value that can be attained by wy , we can show that 


| Mess | < m,Q , k = 0, § = 0, (33) 


m= | ” | | dFyy(a). (34) 


m, 18 called the kth absolute moment of wy . 
If the first K moments are used in estimating F,,(a) from (31), the 
truncation error 7'x is given by 


Te = pe exp (—a'/20") 2 (—)'Hs-(—4/0V2)(1/ov9)'mu/k!. (85) 
T =K+1 

Since it can be shown’ that 
| H(t) | < b2°?-Vn! exp (/2), —-b & 1.086435, (36) 


one can show from (85)-(36) that 


ERROR BOUNDS IN DIGITAL SYSTEMS 3135 


2) MK (Q/c) 25 
o (K+ 1)VK! 
1 = | te 
oVK +14 coVK4+1 
From (31) and (37) one may observe that F,,(a@) may be estimated 
with as great an accuracy as desired if the range of w, is finite and if 
the moments of wy are known. 
If wy is an even random variable, we can also show that 


Mor-1 = O, k21 (38) 





[Tx | < (b/ V/ 2) exp (—a’/4o 


aod; (37) 


10,18 


F,,(a) = 4 erfe (—a/ov2) + vA exp (—a’/20’) 


YY Has-a(—a/oV2)(1/oV2¥*nor/2K)! (89) 





— a2 JA?) BE ——_(@/ oy" _ 
exp (ade) ok OK +2) VOR + 1)! 


[1 — (Q/c)?/ VQK + 2)2K + 3))", 
(2/o)°/[(2K + 2)(2K + 3)? <1. (40) 


| Tox | < 


b 
WV Qe 


By using the inequality” 
| Hoii() | S | t| exp (7/2)(2k + 2)'/(k + 1)! (41) 
we can also show that 


| a | 


oV 





exp (—a"/to?) Ba ALIS LL — (/a)'/(K + 2)", 


(Q/0)°/(K + 2) <1. (42) 
If zz is an even random variable, we have 
F(a) = 4erfe ((—a + toy + Yse)/0V2))on.n 
+ 2ierfe [(—@ + tay — Ysn)/V2))onex- (43) 
Since one can show (see Appendix A) that 


4 erfe (cx + A) + 3 erfe (2 — d) = erfe (2), 5 


[Tox |S 


IV 
oO 


(44) 
we can write 


F(a) = Z(erfe [(—a + Bin) OV 2 ae ’ Sse von = 0, 


= F,,(a), —a+ 2%, 2 0, Yo fae (45) 
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Now from (28) we can write 


F@) = 2(Nexen) onsen 3 (46) 
where 
Nwy.en = erfe (a), 1 = (—a + toy + Yen)/oV2, 
a, exp {—(s + y.,/oV2)"} ds, x2 = (-a+ xy )/ov2. 
oe, (47) 
Since 
exp [—(y.,/oV2)] $1, Vue, (48) 
we have 
ene i: exp [—s° — s(V2/c)yz] ds, (49) 
Var Jes 
F.@) < = ( i * exp [—8 — 8(V2/0) yee] i) (50) 
Via \es Sinus 
7 va (f exp (—s’)&,,(—sV2/c) as) ; (51) 
where 
b= [exp (ty) dF) (52) 


is the moment-generating function of the random variable Zz, . 
If we can find two numbers m, and o? such that 


®, , (t) = exp [tmz + opt /2), V t, (53) 
one can show from (51) that 


F(a) S B.y(@, me , on) = (1 — 02/0)” exp [me/{20°(1 — o2/o°)}] 


el 7a + me/( = op/o) + tus) 27.2 
5 (crf | o2(] a ot/o) ee or/o <1. (54) 





The derivation of the upper bound in (54) is based on results given 
in Ref. 20. 
In this case we then have 


F,,(a) = F(a) = B,,(a, Mr ; oR) 


—a+t+ «ty, 2 0, Gao. A ae ee (55) 
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Since the lower bound in (55) may not be valid if —a + 2, can be 
nonpositive for some value of z,,, , and if the maximum absolute value 
of x,, is 2 monotone-increasing function of N, we note that there is 
an upper bound N,,,.x to N that can be used in estimating the lower 
bound in (55). If this upper bound JN, < ©, we may not be able to 
estimate F(a) from (55) with arbitrarily small error. However if 
there is no finite upper bound to NV such that —a + «,,, 1s nonpositive 
(system with an “open eye pattern’’) and if | mz | and oc, are strictly 
monotone-decreasing functions of N, it is clear that we can estimate 
F(a) from (55) with any desired accuracy. 


III. BOUNDS ON THE TAILS OF PROBABILITY DISTRIBUTIONS 


To use the bounds given in (18), it is necessary to determine L, , (Al) 
and U,,(Au). There are several methods (including numerical methods) 
of determining these parameters, and here we shall discuss two of them. 





From Chebyshev-Bienayme bounds’’’” we have 
Pree S —Al) 5 Gaye, (56) 
Pr fen > Au] s ih, (57) 
where 
(Honor = Yer)» (58) 
Hence we can set 
L,,(—a) = U,, (a) = {Hoplen, (59) 


Also in communication problems, bounds of the Chernoff type have 
been used on the tails of the probability distributions, and these Chernoff 
bounds are often tighter than the Chebyshev-Bienayme bounds.’~°'?*~7° 

One can show” that | 


Pr [2k S —Al] S exp (—) Alexp (—)y- 2) 
= exp (—) Al)&,,(—)), A 2 0, (60) 
Pr [2g > Au] S exp (—A Au)®,,(), = 0. (61) 


The parameter ) is arbitrary and is chosen so as to optimize the bounds 
in (60) and (61). . 
If we can find two functions y,,(—A) and W,,(\) such that 


Os ®, ,(—A) = Wer(—A); Noe 0, (62) 
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and 


0s ,,(\) S ¥.,(A), 2 0, (63) 


and then optimize exp (—A Al)y,,(—A) and exp (—A Au)WV,,(A), we 
can make 


L,,(— Al) = exp (=Yopt Al) W.2(—Nopt) (64) 


and 
U, R (Au) = exp ‘ee opt Au), R (one) : (65) 


The functions y,,(—A) and W,,(A) are often chosen so that (64) and 
(65) have the desired functional forms for optimization.’'?'** From 
(52) one may note that it is not necessary to determine (explicitly) 
®,,(A) to get ¥,,(—A) and W,,(A). Bounds can be used to determine 
these functions. Also one may make use of the semi-invariant moment- 
generating function of zz in determining V,,(—A) and y,,(A). 

If zz is an even random variable, note also that 


®,,(—A) = B2(4), AZ 0, (66) 

and we can make 
en(—) = 2,0), AZO, (67) 
L,,(—a@) = Uzg(a) = exp (— adopt). 2 ont): (68) 


IV. ERROR BOUNDS WITH INTERSYMBOL INTERFERENCE AND ADDITIVE 
GAUSSIAN NOISE 


The methods presented in Section II are now applied to the analysis 
of a binary coherent digital system subject to intersymbol interference 
and additive gaussian noise. Various methods have been proposed to 
evaluate this error probability.’”** They provide either an upper bound 
to the error rate or error rate with a finite pulse train approximation. 

Let us now assume that the signal at the input to the receiver detector 
(see Fig. 3) can be represented as 


Hi) = So aplt — 2) + n(0, (69) 


where n(é) is a gaussian random variable with mean zero and variance o”. 
We shall also assume that {a,} is a sequence of independent random 
variables, and a, = -+1 with equal probability. 

If the zeroth transmitted symbol is a, = 1 and if it is detected by 
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ADDITIVE GAUSSIAN NOISE 







Za, O(t-kT) RECEIVER | 2axp(t-kT) +ni(t) 


FILTER 
R(f) 











CHANNEL 
C(f) 


Fig. 3—Simplified block diagram of a coherent digital communication system. 
Cf), Td), and R(f) denote respectively the transfer functions of the channel, and 


transmitting and receiving filters. 7 is the signaling interval. 


sampling y(t) at t = t) , we can show that 
y(to) = p(to) + De ap(to — IT) + n(to), (70) 


where >,’ does not include the term 7 = 0. Assuming that the slicing 
level of the system is zero, and that there are no other imperfections in 
the system, we can show that the probability of error P, can be written as 


P, = Pr[n + 0’ ap, < —pol; (71) 
where 
Pi = | p(t — IT) |, (72) 
and 
n = n(to). (73) 


Without loss of generality we shall now reorder sequence {p,} in 
such a way that the terms of the sequence are nonincreasing with 
increasing J, and let us denote this new sequence by {r,}. Hence we can 
write 


Gee ee ta E + py arn, < -n| (74) 
or 
P, = Pr[z < —po] = F.(—Po), (75) 
Z2=nt Daan, = ay + ee (76) 
ay =m Do ate, tw = {1, 2,3, --- , N}, (77) 
Wy = 2@y — Nn, (78) 
aa a QT, (79) 


ketne 
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Since zy and Zz, are statistically independent random variables, (18) 
gives bounds to F,(—po). Let us first determine F,,(a) where a = 
—po + Au or a = —pyo — Al. Methods given in Section 2.2 can be 
used in determining F,,(a).* We would like to note here that (31) 
must be used in determining F,,(a) when +1 and —1 do not occur 
with equal probability. 

The recurrence relation given in Ref. 10 to calculate the even order 
moments p2,’S is to be used with care since the summation in the recur- 
rence relation contains both positive and negative terms. In Appendix B 
we give another recurrence relation to compute p,’s (and p2;,,’s, 1 = 0). 
Since the new recurrence relation for po,’s contains the summation of 
positive terms only, we consider this method of computing p,,’s pref- 
erable to that given in Ref. 10. 

We used (31) and our new method for computing p2,’s to calculate 
F,, (a). 

We shall now determine L,,(— Al) and U,,(Au). Since zp is an 
even random variable, we will set Au = Al, L,,(—Au) = U,,(Au). 
Also one can show’ that 


®, (A) I] cosh Nk 


kefn° 


2 
exp I, it + *~ >» “| : A+ A°= ty. (80) 


IA 


From (68) we have 


{Au = >> ri}? 
U,,(Au) = exp | —————=—4— |], Aw Dn Z0,AC fh. 
2 Sr A (81) 
A¢ 
Equation (18) now yields 
[Au —= Sri al 
Fy(—Do — Au)) 1 — 2 exp) — 
[2Ea | 
[Au — >> r,]? 
< FA—p) S F.y(—po + Au) + exp | — : "3 (82) 
: l 


Ae 


For any given N, an optimum Aw can be chosen to minimize the 
difference between the upper and lower bounds in (82). This is often 


* Other methods (including simulation) can also be used in determining F.,(@). 
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found to be difficult and tedious and relies heavily on the search methods 
given in Ref. 7. 
Here we assume that 
A= {x5 (83) 
and we write 


F,,(—po — Au)[1 — 2 exp {—(Au)’/26%}] 
< F(—po) S F.y(—po + Au) + exp [—(Au)?/26n], (84) 
ee Dari. (85) 


Note that any number 62 = ar rv; can be used in computing the 
bounds in (84). This may be done to simplify computing >>, 1°. 

The difference Dy(Au, Au) between the upper and lower bounds 
can be written as 


Dy(Au, Au) = Pr [—po — Au < zy S —po + Au] 
+ exp [—(Au)?/26z]{1 + 2F.,(—po — Au)}. (86) 


Since 6; is a strictly monotone-decreasing function of N, Dy(Au, Au) 
can be made smaller than any given number e. Hence we can calculate 
F',(—po) from (84). 

Several different algorithms can be developed to compute F',(—o). 
One of our algorithms is as follows. Let us assume that we have to 
calculate F,(— po) with a fractional error less than e, . 

Since F,(—po) S 1, we assume that there exists an N such that 


[SE ga Po) =F vay Po) | es (87) 
where 
é S 36 min {F,,(—po), F.y4,(—Po)}- (88) 
For this N we calculate 6; and choose Au so that 
exp [—(Au)*/28z] = 62/3. (89) 
We then calculate Dy(Au, Au) and compare it with 
X, = 4@F.y(—po — Au)[l — 2 exp {—(Au)’/26}]. (90) 


We increase N so that 


Dy. (Au, Au) < xX,! 


N ? 


N'= N. (91) 


It is not necessary to increase N in steps of one. The step size can be 
chosen to suit particular examples. 
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From (18) and (91) we can write 


Ay (—po) S F.(—po) S By (—po), (92) 
Ay (—Dpo) = Fy (—po — Au)[1 — 2 exp {—(Au)’/262}], (93) 
By:(—po) = Fy (—po + Au) + exp {—(Au)’/2B%}, (94) 
By(—po) — An (—Po) S «An (—po). (95) 


It is evident from (92) and (95) that F',(—po) is equal to Ay (— 70) 
or By: (—po) with an error less « . 


We have programmed this algorithm on a digital computer and we 
have been very successful in evaluating F,(—po) from this algorithm. 


4.1 Applications 


Let us now assume that p(t) is obtained by passing a square pulse 
through a single-pole RC-filter or that 


p(t) = 0, t<0 (96) 
p(t) = 1 — exp (—27W)), Ost KT, (97) 
p(t) = exp [—27W(t — T)] — exp [—27rWt], t=T. (98) 


For this pulse we can write 


po = 1 — exp (—27Wt), 0s%4s T, (99) 
and 
r, = [1 — exp (—2rWT)] exp [—27W{t + (k — 1)T}], 
k2=1. (100) 


For 2WT = 0.5, and t = T, we plot in Fig. 4 F,(—po) with an 
error less than 0.2 percent. In this figure we also plot N’, the number 
of terms required in estimating F',(— po). F'.1(—Do) is calculated from 
(31) with a truncation error of less than 0.01 percent. 

Let us now consider the ideal bandlimited pulse p(t) where 





_ sin t/T 
DY) = oe» (101) 
Po — ea ) 6 = oT’ < L, to = 0, (102) 
sin 76 k=, (103) 
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NUMBER OF TERMS N' 


PROBABILITY OF ERROR 
oO. 8 
' i] 
w > 





SIGNAL ~TO-NOISE RATIO IN DECIBELS 


Fig. 4—Probability of error of binary coherent digital system with intersymbol 
interference and additive gaussian noise. The received pulse is an exponential pulse, 
and 2WT = 0.5. The upper bound By:(—po) is plotted in this figure and N was 
increased in steps of one. [By:(—po) —F'.(—po)|/F.(—po) < 0.002. The truncation 
error is less than 0.01 percent. 


sin 76 
ee ee > 
Tar (ke + 6] ; | peas Led (104) 


We shall assume that we take an even number of terms in wy in estimat- 
ing Fy(—Do). 








We have 
oo =. e 
1=2N4+1 
_ sin’ 76 | 1 1 | 
k=N+1 — r (k — 6)” (k + 5)? 
(1 + 8°) sin? 16 E v ‘| 2 
< fiat, iia nee tieee SE pe = 
<2 i= oe le > 1/1 ap (105) 


Since a? is more easily computed than 87, we shall use ag in (84). 

For 6 = 0.05 we plot in Fig. 5, F.(—po) with an error less than 50 
percent when F,(—po) = 2 X 107° and less than 100 percent when 
F(—po) < 2 X 10°°. In this figure we also plot N’ the number of 
terms required in estimating F,(—po). Since a2 is a slowly decreasing 
function of N, the number of terms required for estimating F',(— po) 
is much larger than that in the earlier example. 

Since zy contains a gaussian random variable and since Zz 1S an even 
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random variable, (55) can also be used to obtain upper and lower 
bounds to F,(—po). Equations (53) and (80) can be shown to yield 


Nr = Yr, } (106) 
A 

a = Don: (107) 
Ac 


By choosing A° = ¢,, we obtain the bounds given in Ref. 20. 

Here we would like to note that the relative merits of the two sets 
of bounds cannot be compared as the bounds in (55) may not be appli- 
cable when the system has a closed eye pattern. The lower bound in 
(55) can be shown to be tighter than that in (18) but is not applicable 
to a system with a noneven zz. The random variable zz is noneven 
if +1 and —1 do not occur with equal probability. From the point of 
view of computation, tightness, and applicability, we think that specific 
problems should determine the set of bounds best suited to them. 

The extension of this analysis to m-ary coherent digital systems, 
m > 2, and binary coherent phase-shift keyed systems is obvious from 


NUMBER OF TERMS N'x 1073 
Eo FF FF 72 74 TE IB 8 OR 


PROBABILITY OF ERROR 
cS) 
I 
rs 





SIGNAL-TO-NOISE RATIO IN DECIBELS 


Fig. 5—Probability of error of binary coherent digital system with intersymbol 
interference and additive gaussian noise. The received pulse is an ideal bandlimited 
pulse, and it is sampled at to, fo7’ = 0.05. The upper bound By/(—po) is plotted in 
this figure and N was increased in steps of 100. [By-(—p0)-F'(—po)]/F.(—po) < 
0.5, F.(—po) => 2 X 107%, [By (—po)—F(—po)]l/FA—po) < 1, FA—po) < 2 X 
10~*. The truncation error is less than 0.1 percent. 
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Fig. 6—Distribution function F,(a) = Prizy + zp < al. 
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Refs. 7 and 9. The analysis for higher-order phase-shift keyed systems 
needs extensive modification and will be treated in a future publication. 


V. DISTRIBUTION FUNCTION F',(a) WITH ARBITRARY 2y AND Zp 


Consider two one-dimensional random variables zy and zz. The 
joint probability distribution of zy and 2, is a distribution in ®’, or a 
two-dimensional] distribution. 

Now the probability distribution of 2 = zy + Zp is given by (see 
Fig. 6) 


F(a) = Pr [zy + ze S al 
=[aPS, @weQ if a+ysa (108) 
Q 
and P(S) is the probability function of zy and zz.’ P(S) represents 
the probability of the relation (z, y) C S. 
Since dP(S) = 0, note that (see Fig. 7) 


[ we < is aP(s) + [ APCS) Cea e Oe aoe 


lA 


a+ Al, 


IA 


(z,yeQ. if y (109) 
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[ aP(S) = Prkey Sat All = Fiat Ad, 
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Qe 


F(a) S F.,(a + Al) + Pr [zz S — All. 
Also we have (see Fig. 8) 


(110) 


(111) 


(112) 


/ dP(S) = [ aP(S) ~f PS. Cape. eae 
Q Qs Q4 


Au. 


IV 


(a, y) © Qs if Y 


Since 


dP(S) = Pr [zy S a — Au] = F,,(a — Au), 


Qs 


[ dP(S) 


4 


Pr [er 2 Au], 


F’,(a) 


IV 


F,,(a — Au) — Pr [zg = Au}. 


(113) 


(114) 


(115) 


(116) 
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From (112) and (116) we can write 
F,,(a — Au) — Pr [zn 2 Au] 


< Fa) S F,,(a+ Al) + Prleg S$ —All. (1117) 


Equation (117) is valid even when zy and 2, are statistically dependent 
random variables. 

If the distribution of zp is very much concentrated around some 
point y = Yo, it was shown in Sections II and IV that F,(a) can be 
evaluated with arbitrarily small error if zy and z,z are statistically in- 
dependent random variables and if we can bound F,,(A). If zy and zp 
are statistically dependent random variables, equation (117) shows 
that the same techniques can be used to compute F(a) if the distribu- 
tion of Zp is very much concentrated around some point y = Y%. 


VI. CONCLUSIONS 


We have presented simple upper and lower bounds on the distribution 
function of the sum of two random variables in terms of the marginal 
distribution functions of the variables. 

We have also derived several other bounds when one of the random 
variables is a gaussian random variable or when one of the distribution 
functions is convex or concave. 
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These bounds are then applied to the error rate analysis of a binary 
coherent digital system subject to intersymbol interference and additive 
gaussian noise. Since the difference between the upper and lower 
bounds is a monotone decreasing function of the number of pulses in 
the finite pulse train, the bounds can be used to compute the error 
probability with arbitrarily small error. Application of these bounds 
is illustrated by two examples. Relative merits of the bounds are also 
briefly discussed. 

Many other applications including the analysis of co-channel and 
adjacent channel interference in communication systems will be evident 
to the reader. Some such novel applications will be given in a future 
publication. 


APPENDIX A 
Let us write 
G(a) = 4 erfe (x4 + a) + $ erfe (& — a). (118) 
If x = 0, one can easily show that 
4 erfe (a) + 1 erfe (—a) = 1 = erfe (2). (119) 


We shall now assume that x ¥ 0. Also since G(a) is an even function 
of a, we shall consider a = 0. From (31) and (118) we can write 


1 
G’(@) = —= [exp {—(@ — a)’} — exp {—@+a)"}]. (120) 
Vr 
Note that G’(0) = 0 and that there are no other finite stationary 
points of G(a), « ~ 0. Further one can show that 
G’(a) > 0, x> 0, a> 0. (121) 
G'(a) < 0, x <0, a> 0. (122) 
From (118), (121), and (122) we then have 
4 erfe (cx + a) + f erfe (x — a) 2 erfe (2), 


> x 
4 erfe (ec + a) + 3 erfe (x — a) S erfe (2), e 


IIA. IV 


0, (123) 
0. (124) 
For the sake of completeness we would like to note here that erfe (x) 


is a convex function for x = O and is concave for x < 0. Hence we 
can also show that 


p erfe (x + a) + (1 — p) erfe (@ — a) 
2erfe(r), ttaZz0, t—a20, OSp8il, (125) 
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and 
» erfe (x + a) + (1 — p:p) erfe & — a) 

Serfe(z), t+tas0, t-—-aS0, OSpS1. (126) 
Observe that (125) is not sufficient to prove (123). 


APPENDIX B 


Let 7, denote the partial sum oo £; where 


Wy = 24 g; ; (127) 


and £,’s are statistically independent random variables. From (32) and 
(127) we can write 


Hn = (iy) = 8,(N) (128) 
where 
6,0) = (ni), nel, Oz) = 1. (129) 
Now 
Ok) = (m1 + el"), ok > I, (130) 
or 
alk) = 2 (") a0 —Lonp(h), b> 1, (131) 
where 
ant) = (8), ok) = 1, BEL (132) 
Since 
6,(1) = (mt) = (&) = @,(1), (133) 


and since we shall assume that all a,_,(k)’s are known or can be eval- 
uated, we have a recurrence relation in (131) to compute py, . 

Often &,’s are even random variables, and in this case we can show 
that 


Mot41 = 0, l = 0, (134) 
Mon = Oon(N), (135) 


Ban(K) = >> (5 oat a 1)e2n—2p(K). (136) 


p=0 
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The recurrence relation (136) contains only the sum of positive terms, 
and hence can easily be used to compute uz,’s. 

In Section IV, & = aur: , Mosi = 0,1 = O and a,(k) = rz‘, 7 = 0. 
All even order moments of wy can therefore be easily calculated from 
(136). 

In Refs. 18 and 25 methods have been developed to calculate 2, of 
the random variable @ where 


K 
6 = >. R; cos 6; (137) 
j=l 


and @;’s are independently distributed random variables uniformly 
distributed over the range [0, 27). Most of these methods use an infinite 
series expansion, and often the accuracy obtained from these methods 
is questionable.” 

Noting that we can set £; = R; cos 0; , Mor41 = 0,1 = O, and 


a2i(j) = ((R; cos 6;)**) 
or 


(26)! 
2”*(¢!)? ] 





ai(j) = Rj (138) 


all even order moments y2,’s can be calculated by using (136) and (188). 
This method of calculating u2,’s can be shown to be analogous to 
that given in Ref. 26 and is preferable to that in Refs. 18 and 25. 
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Optical Fiber J oining Technique 


By D. L. BISBEE 
(Manuscript received July 15, 1971) 


This paper describes a method of thermally fusing clad glass fibers, end 
to end, to obtain a good mechanical joint with low transmission loss. 
Methods of preparing fiber ends and aligning them for joining are dis- 
cussed. Two sizes of fibers were joined (10.8-ym core and 20-um core clad 
fibers with outside diameters of 75 wm and 150 wm respectively) .* 


I. INTRODUCTION 


There is a great deal of interest in using glass fibers as optical wave- 
guides to carry information in much the same way as wires or metallic 
waveguides do. If glass fibers are to be used in this way, they will 
need to be joined just as wires and metallic waveguides must be joined. 

A method for joining single fibers was developed. Clad glass fibers 
were joined which had cores of 10.8-um and 20-um diameters and 
overall diameters of 75 um and 150 um respectively. The cores were 
Schott SSK-1 glass and the cladding of Schott SK-14 glass which 
have glass transition temperatures of 621°C and 649°C respectively. 
Good mechanical joints which can be made quickly with transmission 
losses as low as 11.5 percent. were obtained, but lower losses should 
be possible with a little more effort. 


II, FIBER END PREPARATION 


To get a good joint, good fiber ends are needed. Polishing or etching 
the fiber ends has been suggested, but we have found that if a fiber 
is broken properly it will have an end that is suitably flat over most 
of its surface and perpendicular to the axis of the fiber as seen under a 
microscope. Figure 1 shows two good ends of 10.8-um core, 75-ym o.d. 
fibers magnified 500X. The break can be made by scoring the fiber 
with a razor blade and breaking it or by laying the fiber across a sharp 
metallic edge and positioning a Tesla coil so that its discharge is con- 


* The fibers were manufactured by DeBell and Richardson, Inc., of Hazardville, 
Connecticut. 
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Fig. 1—Two cleanly broken ends of 10.8-um core, 75-um o.d. fibers magnified 500X. 


centrated at the point where the fiber touches the metal, then breaking 
it by bending or pulling. The sparking must continue for several seconds 
with greater time for greater fiber diameter. To prevent the deposition 
of metal on the fiber because of the electrical discharge, platinum wire 
was used as the metallic contact. Fibers of 250-um o.d. or greater can 
be scored with a file or diamond before breaking. 

If the fiber is bent to the breaking point after it 1s scored or treated 
with an electrical discharge, one end will have a lip protruding and the 
other end will have a corresponding absence of material. This can be 
seen in Fig. 2. This lip is sometimes as long as the diameter of the 
fiber and would prevent the end from being brought close enough to 
another fiber end to permit joming them. If the fiber is pulled instead 
of bent until it breaks, the lip is not produced on most occasions, so 
this is the recommended procedure. 

Very small fibers will sometimes soften and bend from the heat 
when the discharge from the Tesla coil strikes them, but with care, 
fibers as small as 25-um o.d. have been broken with good ends suitable 
for joining. 
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WI. ALIGNMENT 


After the ends of the fiber have been prepared, they must be mounted 
so that they can be aligned. Teflon-coated tweezers to hold the fiber 
ends were mounted, one on a general purpose 3-dimensional micro- 
manipulator and one on a precision 3-dimensional micromanipulator 
with a positioning resolution of 0.127 um. 

The Teflon-coated tweezers are gentle with the fiber and allow it to 
slip when the fiber contracts after having been heated, as will be dis- 
cussed later. 

From measurements of light output versus fiber end displacement 
we find that a misalignment of less than 2 ym in the 10.8-ym core fiber 
gives 10 percent less transmission than when the ends are aligned. 
Losses due to fiber offsets are covered in detail in Ref. 1. 

One can determine when the ends are aligned by viewing them through 
a microscope and assuming that the core is concentric with the outside 
of the cladding. To get two perpendicular views of the fiber, a mirror 
can. be mounted so the fiber can be viewed directly from the front and 





Fig. 2—A characteristic break in a 10.8-um core, 75-um o.d. fiber showing the 
protruding lip. 
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through the mirror from the side. Two problems in this method of 
alignment are, first, the core and cladding may not be concentric, 
and second, rather sophisticated optics are needed to see an alignment 
error of 1 wm or less. If the core and cladding are made concentric and 
if a microscope of 200X or greater is used, one can probably align them 
well enough by this means. 

Another method we have used of determining optimum alignment 
is to send laser light down the fiber to a detector and adjust the fiber 
ends for maximum transmission. A problem with this method is that 
if the end of the fiber is broken at an angle with respect to the normal to 
the fiber axis, the maximum transmission will be obtained when the 
fiber ends are misaligned to compensate for the offset in the beam 
direction caused by the angle of refraction at the nonnormal surface. 
This offset is y in Fig. 3. This error is small, though, if the angle is 
small. With a fiber of 80-um diameter, core index of 1.6, and a surface 
at the end that is 10 degrees from the normal, the ends of the fibers would 
have to be offset 0.73 um to correct for the beam misalignment. This 
would introduce about 5 percent loss in our 10.8-um core fiber. An 
angle of 10 degrees is large, so one should be able to do much better 
than that. 


IV. JOINING THE FIBERS 


Several unsuccessful attempts were made to join fibers with epoxy 
alone and epoxy in a glass sleeve. When using epoxy alone, the resultant 
joint was too weak to keep the fiber ends aligned when transverse 
pressure was applied. When using sleeves, the tolerance between sleeve 
id. and fiber o.d. had to be very close, of the order of 1 um to keep 
the fibers aligned properly, and this tolerance is hard to obtain. Further, 
a bubble formed at the fiber junction inside several sleeve joints. Thus, 
this method was considered limited in practicality. 

A method that worked was the fusing of the fiber ends. Number 
24 nichrome wire was wound around two metal posts so as to leave an 
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Fig. 3—Beam refraction at a nonnormal fiber end. 
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Fig. 4—Heating wire with fiber ends held in place by Teflon coated tweezers. 


opening between the wires through which the fiber ends could pass 
as shown in Fig. 4. 

The fiber ends were aligned, leaving a space of about 15 um between 
them for thermal expansion. About 14 amperes of current were passed 
through the nichrome wire which surrounded the fiber ends causing 
the wire to heat up and fuse the fiber ends together. The longitudinal 
expansion of the fiber when heated closed the 15-um gap that was 
left between the fiber ends. Of course, when the fiber cooled it shrank 
again, but the fiber could slip in the Teflon-coated tweezers when 
shrinking. The ends fused together in about 30 seconds after the heat 
was applied. To tell when the ends were fused, a lamp was placed so 
that the specular reflection from the fiber ends could be seen in the 
microscope. Disappearance of the reflection indicated the surface had 
vanished and the ends were fused. Figure 5 is a microphotograph of a 
fused joint in a 10.8-um core, 75-um o.d. fiber, at 500X magnification. 
Such joints exhibited losses as low as 11.5 percent. 
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JOINT 





Fig. 5—A thermally fused joint in 10.8-um core, 75-um o.d. fiber. 


V. CONCLUSIONS 


After two fiber ends are carefully broken and properly aligned, they 
can be fused together by means of a heated wire to give a good me- 
chanical joint with an acceptable amount of loss. Fibers of 10.8-um 
core were joined and gave losses as low as 11.5 percent, but with a 
little more effort it is believed that joints could be made with con- 
siderably lower loss than this. 
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Measurements of Loss Due to Offsets 


and End Separations of Optical Fibers 


By D. L. BISBEE 
(Manuscript received June 23, 1971) 


If fibers are to be coupled together by means of detachable connectors, 
there 1s a need to know how much light will be lost by misalignment or 
axial separation of the fiber ends. 

Measurements were made of coupling efficiency from one fiber to another 
versus offset and end separation with and without tindex-matching liquid 
between the ends for a single-mode and a multimode fiber at X = 0.6328 ym. 
Graphs are presented for offsets as great as 3 radii and for end separations 
up to 127 um. Maximum coupling efficiency of 97 percent was obtained, 
and about 50 percent was obtained with an offset of 1 radius. 


I. INTRODUCTION 


Coupling of glass fibers, end to end, for transmitting optical signals 
requires critical alignment at each coupling point. Coupling can be 
accomplished by permanently fusing the fibers together,’ by using 
permanently bonded sleeves, or by using detachable connectors. 

This paper presents experimental data for the coupling efficiency of 
light from one fiber end into another as a function of offsets and axial 
separation which will be present in detachable connectors and to a 
lesser extent in fused joints and sleeve joints. 

Measurements of coupling efficiency versus end separation and 
translation were made and plotted for a multimode fiber with core 
diameter of 10.8 um and a single-mode fiber with a core diameter of 
3.7 pm.* 

Several people’'*"* have studied the problem of launching efficiency 
from a Gaussian beam into a fiber, both theoretically and experimentally. 
A comparison will be made between the theory and our experimental 
data. 

* The 10.8-um core fiber was manufactured by DeBell and Richardson, Inc., 


a Hazardville, Connecticut, and the 3.7-um core fiber was made by Corning Glass 
0. 
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Il. EXPERIMENT 


2.1 Hqutpment 


Figure 1 is a sketch of the measuring apparatus. The light from a 
He-Ne laser operating in the single transverse TEMo) mode that gives 
a Gaussian distribution was chopped, then focused down by a 10X 
microscope objective into a fiber that we call the launching fiber. The 
output end of the launching fiber was mounted on a precision 3-dimen- 
sional micromanipulator with positioning resolution of 0.127 um. The 
output of this fiber was launched into a fiber that we call the receiving 
fiber. The output of the receiving fiber was detected by a solar cell 
immersed in index-matching liquid, measured by a lock-in amplifier, 
and recorded on a chart recorder. Mode strippers (not shown) were 
used on both launching and receiving fibers to eliminate the light 
traveling in the cladding. These were made by cutting an S-shaped 
groove about 0.25 mm wide and deep in a piece of plexiglass about 
30 cm long and placing the fiber in the groove, then filling it with index- 
matching oil of a slightly higher index than the fiber cladding. None of 
these dimensions are critical, and the same effect could probably be 
obtained by making the stripper much smaller. 


2.2 Method 


(z) Coupling efficiency versus fiber end separation was measured by 
first aligning the fiber ends axially and then bringing them into contact 


O 
B Fe 


Fig. 1—Schematic diagram of light coupling measuring apparatus: B, chopped 
beam from laser; O, microscope objective; F;, launching fiber; F,, receiving fiber; 
M, index-matching oil; D, detector. 
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with each other for the measurement at zero separation. One fiber was 
then separated from the other by known amounts and the output of 
the receiving fiber was recorded. This was repeated several times for 
each case and the mean of the readings was plotted. These measure- 
ments were made for the 10.8-um and 3.7-um-diameter core fibers with 
and without index-matching oil in the gap. A reliable measurement 
could not be made at zero separation without oil because of interference 
effects. . 

(iz) Coupling efficiency versus fiber lateral displacement was measured 
and plotted for both fibers with and without matching oil between them. 
These measurements were made for fiber end separations of 5.08 um, 
10.16 pm, 25.4 pm, and 50.8 um. The smallest separation for which 
scanning measurements were made was 5.08 um because they could not 
be made at zero separation as one fiber would rub against the other, and 
also because of the interference effects mentioned above. In each case 
the scanning was begun with the fibers out of alignment on one side of 
the fiber axis so as to give approximately zero transmission, then one 
fiber was translated past the axis of the other reaching a maximum 
transmission, and continuing on until approximately zero transmission 
was reached again. This operation was repeated four times, once each 
in the plus and minus z and y directions of Fig. 1 for each separation. 
The mean of these four readings was then plotted. 

(it) The coupling efficiency expressed in percent was calculated from 
the power levels measured at the ends of the launching and receiving 
fibers and from the transmission loss in the second fiber. Absolute values 
of the data presented could be in error by as much as 3 percent. 


2.3 Results 


2.3.1 Multimode 10.8-um-Diameter Core Fiber 


For this fiber r = 5.4 um, \ = 0.6328 um, the core index is 1.6171, 
and the cladding index is 1.6038. From this data we find that the fiber 
is capable of propagating 67 modes.” ; 

It would be difficult to mathematically predict the light coupling 
loss due to offsets or fiber end separation from such a multimode fiber, 
So no comparison was made between theory and experimental data on 
this fiber. 

The coupling efficiency at 2.54-um separation without oil in the gap 
was measured to be 88.8 percent, and at zero separation with oil in the 
gap it was 97.06 percent. These values appear as points on the curves 
of Fig. 2. The data in graph B was taken with index-matching oil 
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Fig. 2—Curves of coupling efficiency versus fiber end separation for 10.8-1m 
Btgcigiad core fiber: A, without index-matching oil in the gap; B, with matching oil 
in the gap. 


(n = 1.6204) in the gap between the fiber ends, and in graph A without 
oil. A reliable reading could not be taken at zero separation without oil 
in the gap because of the interference between the transmitted and 
reflected beams between the plane parallel surfaces of the fiber ends, 
but all the readings for each graph at 2.54-um (approximately 0.5 
radius) separation were within -+0.5 percent of each other. 

The Fresnel reflection loss at the two core-air interfaces is 11.2 percent 
(neglecting resonance effects). If we subtract this from 100 percent we 
get 88.8 percent which is the coupling efficiency that was measured at 
2.54-um separation without oil in the gap. This implies that all the light 
except that reflected is coupled into the receiving fiber. 

When oil of index n = 1.6204 was put between the fiber ends elimi- 
nating the core-air interfaces, the amount of light coupled increased by 
about 8 percent to 97 percent, which was less than the 11.2 percent 
increase predicted. 

Figure 3 shows the coupling efficiency versus offset for four different 
fiber end separations without index-matching liquid in the gap. Figure 4 
has the same set of curves with index-matching oil in the gap. 

We see that only 50 percent of the light is transmitted with perfect 
axial alignment for fibers separated by 50.8 um without matching oil 
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Fig. 3—Coupling efficiency versus displacement of 10.8-um-diameter core fiber 
without index-matching oil in the gap. 


between them. It can be seen also that all the curves are coincident at 
about 20 percent efficiency at a little less than 1.5 radii translation both 
with and without index-matching liquid in the gap. 

2.3.2 Single-M ode 3.7-um-Diameter Core Fiber 


For this fiber, the fiber characteristic term® 
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Fig. 4—Coupling efficiency versus displacement of 10.8-um-diameter core fiber 
with index-matching oil in the gap. 
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is 2.2 at \ = 0.6328 um, where ais the radius of the core, 7, is the refrac- 
tive index of the core, n, is the index of the cladding, and \ is the wave- 
length of the light. Only the HE,,; mode can propagate because all other 
modes are cut off in a fiber with R < 2.405. 

Coupling efficiency versus fiber end separation is plotted in Fig. 5. 
The coupling efficiency at 2.54-um separation without oil in the gap was 
found to be 90 percent. This is shown in curve A, and when this curve 
is extrapolated to zero separation, we find the coupling efficiency at 
zero separation is approximately 91 percent. This compares with 90 
percent for the multimode fiber at zero separation without oil. With oil 
in the gap at zero separation the coupling efficiency was 97 percent. 
Thus, the matching oil increased efficiency by 6 percent. This is the same 
efficiency as was found for the multimode fiber with oil in the gap at 
zero separation. 

There are no striking differences between the coupling efficiency 
versus fiber end separation curves of the two fibers which are shown 
in Figs. 2 and 5. It has been noted that at zero separation the efficiencies 
are virtually the same, but at 127-um separation, which is almost 
70 radii for the 3.7-um-diameter core fiber though less than 24 radii 
for the 10.8-um-diameter core fiber, we see some small differences. 
Without oil at this separation the 10.8-um core fiber has 14.5 percent 
efficiency and the 3.7-um core fiber has 13 percent showing the multi- 
mode fiber more efficient by 1.5 percent. With index-matching oil 
at this separation, coupling into the single-mode fiber is more efficient 
by 5 percent, than into the multimode fiber. At about 50-um separation 
both with and without oil in the gap, coupling efficiency is about 10 per- 
cent less in the single-mode fiber than in the multimode fiber. 

Figures 6 and 7 show curves of coupling efficiency versus fiber trans- 
lation for the single-mode fiber. Figure 6 shows curves of coupling 
efficiency versus fiber translation for four end separations without 
matching oil in the gap and Fig. 7 shows the same curves with oil in 
the gap. These curves are comparable to the corresponding ones for the 
multimode fiber in Figs. 3 and 4, except for the above-noted lower 
efficiency for the single-mode fiber at greater end separations. We see 
also that the efficiency does not drop to zero as quickly with respect to 
displacement with the single-mode fiber as with the multimode one. 
This is understandable in that with the single-mode fiber, the field 
extends into the cladding, but in the other fiber the field is rather sharply 
restricted to the core. 
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Fig. 5—Coupling efficiency versus fiber end separation of 3.7-um-diameter core 
fiber: A, without index-matching oil in the gap; B, with index-matching oil. 
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Fig. 6—Coupling efficiency versus displacement of 3.7-~m-diameter core fiber 
without index-matching oil in the gap. 


3166 THE BELL SYSTEM TECHNICAL JOURNAL, DECEMBER 1971 


III. COMPARISON WITH THEORY 


In the absence of irregularities such as nonflat fiber ends and index 
mismatch, the coupling efficiency from one fiber to another at zero 
separation should be 100 percent. 

Launching efficiencies of nearly 100 percent for the HE,, mode from 
an incident Gaussian beam have been predicted theoretically by J. R. 
Stern and R. B. Dyott,’ and D. Marcuse.* The field profile of the HE,, 
mode is similar to a Gaussian distribution, thus the theory which predicts 
the launching efficiency from a Gaussian beam into a fiber should ap- 
proximately predict the coupling efficiency from our single HE,, mode 
launching fiber into the corresponding receiving fiber. Marcuse cal- 
culated the theoretical curve of coupling efficiency versus translation in 
Fig. 8 using the parameters of our single-mode fiber. The experimental 
data for our fiber with index-matching oil in the gap at 5.08-um separa- 
tion was modified to represent a curve of zero separation and is plotted in 
Fig. 8 for comparison. At the top of the curve we launched only 97 per- 
cent of the incident light where theoretically almost 100 percent could 
have been launched. In the middle the curves are coincident and at the 
lower end the discrepancy between the curves is probably due to the 
fact that the edge of the HE,, mode is not well approximated by the 
Gaussian beam. 
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Fig. 7—Coupling efficiency versus displacement of 3.7-um-diameter core fiber 
with index-matching oil in the gap. 
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Fig. 8—Coupling efficiency versus translation of 3.7-um-diameter core fiber: 
solid line is theoretical, broken line experimental. 


IV. CONCLUSIONS 


With matching oil in the gap between two fiber ends with zero 
separation, 97 percent of the light from one fiber was coupled into the 
other in either single-mode or multimode fibers. The light coupled 
from a single-mode fiber into another single-mode fiber follows closely 
the theory for light coupled from a Gaussian beam into a single-mode 
fiber when one fiber is offset with respect to the other. For fiber ends 
separated up to 25 um, a lateral displacement of one radius reduces the 
coupling efficiency to about 50 percent if there is index-matching liquid 
between the fiber ends, and to about 40 percent if there is none. The 
coupling efficiency was reduced by 10 percent with end separation of 
about 3 radii for the multimode fiber and about 7 radii for the single- 
mode fiber. 
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Zero Loss Transfer Across Gaps in a CCD 


By R. H. KRAMBECK 
(Manuscript received July 19, 1971) 


When a charge-coupled device is made with a single layer of metalliza- 
tion, adjacent electrodes must be placed several micrometers apart. As a 
result there may be some difficulty in moving the charge from one electrode 
to another. In this paper it will be shown that, for any substrate material, 
there 1s a wide range of interface charge for which complete transfer can 
be achieved regardless of electrode separation. It will further be shown 
that for a p-substrate with a doping of less than 10'°/cm’*, the natural 
charge found in a good quality thermally grown layer of SiO, 1s always 
of the appropriate sign and magnitude to ensure complete transfer. There- 
fore, for simple fabrication of a CCD with one layer of metal, this substrate 
material ts the appropriate choice. 


I. INTRODUCTION 


In a charge-coupled device,’ charge is transferred from a potential 
well under one metallization to a deeper potential well under an adjacent 
metal. Since’ the metals cannot touch one another, this process involves 
passing through an interelectrode gap. If the potential pattern is such 
that there is a potential barrier or a potential well in the space at the 
end of transfer, then complete charge transfer is impossible. If two 
layers of metallization are used (which are separated by an insulating 
layer of 1000 A), no barrier or well can form. This requires a more com- 
plex technology, so it would be desirable to make a CCD with one 
layer of metallization. 

In this paper the effect of interface charge and substrate doping on 
potential well and barrier formation will be analyzed. It will be shown 
that, for any given substrate doping and drive voltages, there is an 
interface charge for which neither well nor barrier forms. The calcula- 
tions in this paper will assume an infinitely long gap which will ensure 
that no yield losses will result from small variations in gap width or 
surface charge magnitude. 
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II. ANALYSIS 


The objective of this analysis will be to find the conditions which 
will permit complete charge transfer. It is therefore necessary to examine 
the state of the CCD at the end of the transfer of only the largest ONE 
that can be accommodated. That is, if the largest ONE is transferred 
without loss, then any piece of information could be transferred with 
no loss of charge. Our objective then is to find free-carrier density vs 
position at the end of transfer. If there is neither a peak nor a minimum 
of free-charge density in the space, then the space is not interfering 
with transfer. 

To find free-carrier density on the surface, we must find total surface 
charge density and subtract fixed charge density. The first part of the 
problem is to find charge density vs position at the end of transfer. 
This is simplified by the fact that at the end of transfer of the largest 
possible ONE, the semiconductor surface is an equipotential (the 
movement of free charge will not stop until an equipotential is achieved). 

The geometry of the interelectrode space is shown in Fig. 1. Just 
below the semiconductor surface the field is uniform because the 
surface is an equipotential. The displacement field is given by 


D, = V2Ne.Vq (1) 


where N is the semiconductor doping, e, is the semiconductor dielectric 
constant, V is the surface potential, and q is the electronic charge. 
Just above the surface the field is 


Dy, a D, =f qn — Q.. (2) 


where n is the free-electron density (em ”) and Q,, is the density of 
charged states either in the oxide or at the semiconductor surface 
(p-type substrate has been assumed with electrons as free carriers). 
From equations (1) and (2) 


qn = Dux ey 2Ne,Vq se Ques . (3) 


Since Q,, and V are constants, differences in electron density from 
one point to another are directly related to changes in D,, by 
AD, : 

qd 


An = (4) 

This means that if D,, varies monotonically from the region under 
one electrode to the region under the adjacent electrode, then electron 
density also varies monotonically and no barrier or well exists. In Fig. 2 
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Fig. 1—Cross section of a charge-coupled device showing the interelectrode gap 


the three possible variations of D,, with position are shown. First, in 
Fig. 2a V < V, and V < V2 where V, and V; are the voltages applied 
to the two metals. In this case, D,, is positive under both metals but 
drops to zero in between (for sufficiently large separation). This means 
electron density is lower in the space between the electrodes than it 
is under either of them, and as a result a barrier must be holding back 
charge. Second, in Fig. 2b V > V, and V > V.. Here D,, is negative 
under each electrode but is again zero in the gap. This means electron 
density is higher in the space than under either electrode, so a potential 
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Fig. 2—The variation of oxide field at the surface of the semiconductor, Doz, vs 
position along the surface when the interelectrode gaps are arbitrarily large. The 
curves are for three combinations of applied voltages, Vi and V2, and surface po- 
tential, V: (a) V < Vi, V < Va, (b)V >UN1,V >V2,(c)V 2V1,V S Va. 
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well must be present. Finally, in Fig. 2e V = V, and V S V;. Here 
D.; S$ 0 under metal 1, D,. 2 0 under metal 2, and D,, = 0 some- 
where in between. This is the situation desired. It is therefore necessary 
to find only what value of Q,, will cause the surface potential at the end 
of transfer to be in between the two metal voltages. 

Let us assume metal 1 is giving its charge to metal 2. Then at the 
end of transfer free-electron density under metal 1 is zero. Therefore 
from equation (3) 


Diz, = V2NeVq — Qu. (5) 
Also the difference between V, and V is given by 
V,-V= 25 6) 


where 6 is the oxide thickness. Since we need V = V,, D,., < 0. There- 
fore: 


Q..2 V2Ne,Vq. (7) 
The minimum value of Q,, is obtained when V = V,(D,., = 0) which 
gives: 

Q..2 V2Ne,Vig. (8) 


As an example, with N = 10'*/cm’, e, = 107” F/em, and V, = 1 volt, 
Q,./q = 3.5X10'°/em’. 

This specifies the minimum for N,,. Any lesser value would give 
rise to a barrier. To prevent formation of a well, we have V S V2. Under 
electrode 2 


Diz, = Ds — Qss + Qn 
Diez 2 9. 
Therefore 
Qu < V2Ne.Vq + qn. (9) 


If the largest allowable Q,, is being used, D,., = Oand V = V,. Also, 
from equation (4) 


Doz, ~ Dos, = Qn. 
Therefore 


V¥,.-—V,=8%3 


Cox 
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where e,, 18 the oxide dielectric constant. The maximum value of 
Q,; 18 therefore 


S V2NeVig + “2 (V2 — Vi). (10) 
Using V. = 10 volts and the same numbers as above for the other 
parameters gives 
Qs: 12 2 
; < 2.0 X 10°/cm’. 


Equations (8) and (10) specify the allowable range for Q,, which is 
the result needed. 


III. DISCUSSION 


The previous section showed that a specific range of Q,, is appropriate 
for any given substrate doping and pair of driving voltages. In this 
section the attainability of this range will be discussed. The range is 
shown graphically in Fig. 3 for the driving voltages V, = 1 volt and 
V2 = 10 volts. 

It turns out that for a thermally grown oxide on silicon Q,,/q in 
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g. 3—Minimum and maximum surface charge vs substrate doping for V; = 1 
ale ae for V. = 5 volts and 10 volts. 
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the oxide is normally in the range 10** to 6 X 10°'/cm’. Therefore 
for a p-type substrate any doping under 10'°/cm’* will ensure complete 
transfer. For other substrates and oxides some doping must be in- 
troduced. Probably the easiest way would be a uniform ion implant 
chosen so the total charge including oxide charge lies in the allowed 
range. 

An important point about Fig. 3 is the width of the allowed range. 
For drive voltages of 1 volt and 10 volts the allowed range of charge 
extends over at least one order of magnitude for any substrate doping 
below 10'°/cm® which ensures that a proper amount of charge can be 
obtained even with substantial slice-to-slice variations in charge. It 
also means that aging effects are unlikely to change the amount of 
charge to a value outside of the desired range. 

In the analysis it was assumed there is sufficient time for complete 
transfer to occur. An appropriate question is how much longer transfer 
takes if a gap separates the metals. R. J. Strain and N. L. Schryer’” have 
shown that in a CCD the rate of carrier flow is inversely proportional to 
the sum of the capacitances between surface and substrate and between 
surface and metal. Since this sum is much higher when the surface is. 
covered by a metal, the low capacitance gap will not significantly 
increase transfer time as long as the gap is shorter than the metal. 


IV. CONCLUSIONS 


From the above discussion it can be seen that the simplest way to 
obtain complete transfer in a charge-coupled device, with a thermal 
SiO, layer used as insulator, is to use as substrate material p-type 
silicon with a doping of less than 10°°/em’. 

It can be further concluded that for other substrate materials and 
for other insulators an appropriate charge range always exists, and 
that this charge range can be obtained by ion implantation over the 
entire surface of the slice without use of a photographic masking step. 
Selective implantation of the region between the metals and the regions 
under the edges of the metals could also eliminate potential barriers 
and wells. However, because of fringing fields, this can be achieved 
only by extension of the tail of the charge distribution for a specific 
distance under each metal, but this would be a more difficult way to 
achieve complete transfer. 

Summarizing, it should be clear that potential barriers and wells 
in the spaces between electrodes in charge-coupled devices can be 
eliminated easily. Therefore the possibility of incomplete charge transfer 
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should not be cited as a reason for spacing electrodes abnormally 
closely or using two layers of metallization. For most applications 
future CCD’s should be made on p-type substrates with less than 
10°°/em* doping since the proper charge is automatically obtained. 
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A General Class of Rate-Change Circuits 


By 8S. V. AHAMED 
(Manuscript received June 3, 1971) 


Rate changing ts encountered in data transmission, distributing, collating, 
encoding, and decoding. It entails the transformation of data at one rate 
to data at another predefined rate. It 1s the object of this paper to investigate 
a special class of circuits that accomplish such a transformation quite 
mechanically and methodically when only one clock is available to propa- 
gate all the binary bits of information. These circuits can be implemented 
by magnetic domains, by charge transfer and charge-coupled devices, or by 
any other technology which permits (i) propagation of binary bits of 
information by a certain modular “period” or ‘distance’ in a modular 
unit of time and (2) gating of preselected bits of information from one 
branch into one or the other branches of three branch nodes in a circutt. 
The generalization of rate-changing circuits to distributing, collating, 
reversing, and reversing-with-rate-changing is also presented in the paper. 


I. INTRODUCTION 


Algebraic coding schemes’ increase the length of the information 
block. Shortened codes’’*® obtained from such coding schemes further 
alter the length. When these codes are transmitted uniformly over the 
telephone networks, a change of rate from the received rate to the 
transmitted rate is essential. Further, channels carry coded information 
at different rates. If it is necessary to obtain information at a uniform 
rate, encode it with a general coding scheme and transmit it over any 
general channel, decode it to its original rate and refurnish it, then many 
stages of rate changing are necessary. 

To collate bursts of information from many sources onto a signal 
channel also entails a rate-change process. Distributing circuits which 
select and uniformly distribute a given number of information bits from 
longer blocks of uniformly received information are special cases of 
rate-change circuits. It is the object of this paper to report a set of 
circuits which perform (2) general rate changing, (27) general collating 
and sorting, and finally (277) reversing and reversing-with-rate-changing 
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functions. They operate methodically, and mechanically when only one 
clock source ts available for propagation of information bits. This condition 
exists specifically in magnetic domain technology*”’ with field access 
propagation and to a less stringent measure in charge-transfer’ and 
charge-coupled’* device technologies. In general, the circuits presented 
can be implemented by any technology that permits (z) all binary bits 
of information to be shifted from one location to the adjoining location 
within one clock cycle and (22) the channeling of information bits to one 
or the other branch of a node within a circuit. 

For clarity of exposition, the following terms are defined. A “‘period”’ 
is a unit of distance or location by which a binary bit of information 
may be moved in one unit of time. A “clock cycle” is the smallest unit 
of time. It is also the time required to propagate one binary bit of 
information by one period. A ‘‘gate” is a device which channelizes a 
selected number of information bits into one or another branch at a 
node in a circuit. Such functions can be readily accomplished*’®’’’ in the 
various technologies*’* presently under development in the Bell 
System. | 

The paper is divided into five sections. After the introduction, Section 
II describes a general class of rate-change circuits. These circuits become 
an important element in the collating and distributing circuits described 
in Section III, and Section IV discusses reversing circuits by which the 
order of incoming data bits may be reversed. Section V is devoted to 
conclusions. 


II. GENERAL RATE-CHANGE CIRCUITS 


The principle of rate changing is explained by two examples in this 
section. The generality of the principle is proved in Appendix A. 


2.1 Rate-Reducing Circuits 


Example 1. Consider a block of data 8 bits long to be expanded by 
a ratio of 2, yielding a new rate that is half the original rate. Figure 1 
represents a circuit in which the data arrives at the clock rate. The 
incoming polynomial may be written as | 


Up = Ao + a,X + AX” + a3X° + Q4X* + asX° + gx” -- aX’. 


The four data positions corresponding to the first four terms of uy are 
diverted into the upper part of stage 1 containing N, periods, and the 
second four data positions are diverted in the lower half of stage 1 
containing N, + 4 periods. Gate g, diverts the first two data positions 
into the upper half of stage 2, the next two data positions into the lower 
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BINARY POSITIONS ot 2 1 
DIVERTED ——~ t { { 


N, No Nz, 


SOURCE 





—~<- 


BINARY POSITIONS +4 2 
DIVERTED — STAGE 1 STAGE 2 STAGE 3 


NOTATION: [s>- SOURCE, <s]-sink, (g- GATE 


Fig. 1—A simple example of a 2:1 rate-reducing circuit for an 8-bit data block. 


half, and so on, so that 
uy = X™(ao + aX + aX? + a3X°) 
ee XN *4(a,X* + as,X° + agX° + aX’) 
= X™(xX%q, 4+ a,X + aX? + a,X°) 


+ X*(a, + asX + aX? + a,X°)}. 
Similarly 


Up = X™ 4X" 4X, + aX) + X%(a.X? + aX’) 
+ X°(a,X*° + asX°) + X°(a,X"° + a,X")} 
= X™™ 1 X%ag + aX) + X*(a, + a;X) 
+ X*(a, + a5X) + X (ag + a;X)} 


ae 2a? 2et 2! 20<-—— BINARY 
4 POSITIONS 
i f t DIVERTED 
= Ne Na-1 Na 
aie cs 
BS = U _ 
Uo uy \ rain TAs, a 
\----»/ 
Supe iie: Naver Nait2k Na+k 
! | 
i _— BINARY 
<7 ies as 2" 20“ BOSITIONS 
a= Locs’; k=m-1 DIVERTED 
SOURCE STAGE! STAGE 2 STAGE L STAGE a-1 STAGE a SINK 


Fig. 1a—General rate-reducing circuit for an n-bit data block and rate ratio of 
mil. 
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and 


Us; = tN ats 


{Qo + a,X? + a,X* + a3X° + a,X* + aX? + a,X” + a,X"}. 


The rate of the polynomial u; corresponds to half the rate of the input 
polynomial uy . The output polynomial u; has a delay of (V, + N, + N3) 
clock cycles corresponding to the sum of periods in the upper section 
of the rate-change circuit. In charge-transfer and charge-coupled’® 
devices (VN, + N. + N3;) can be made zero. In magnetic domain circuits 
a certain minimum value for (V, + N. + N;) is foreseen. 

In general it can be seen from Appendix A that the number of stages 
in a rate-reducing circuit is determined by the number of bit positions 
in the data block. The differences in periods between the upper half 
and the lower half of each of the stages are determined by the rate-change 
ratio. Table I contains the characteristics of a general circuit (Fig. 1a) 
for an n-bit data stream. The rate-change ratio is m. In Table I, a, 
denoting the number of stages, is an integer which satisfies the equation 


2°2n or a2 logs, (1) 
and k is an integer obtained as 
k=m-— 1. (2) 


The delay in the circuit can be computed as (Vi + N.+ N;--- N,) 
clock cycles (see Fig. la) from the normal polynomial calculations 
presented in Appendix A. 

2.2 Rate-Increasing Circuits 


Rate-increasing circuits are derived as inversions of rate-reducing 
circuits. If the data flow is reversed in a rate-reduction circuit, a rate 
increase results. The process is investigated in the following example. 


Example 2. A block of data 16 bits long, each bit arriving every four 
clock cycles, is to be condensed to a data block of 16 clock cycles thus 


TABLE I—CHARACTERISTICS OF A GENERAL RATE-REDUCING CriRcUIT 


Ratio of Number of Bit Positions Difference in Periods 











Input Rate Diverted by Gates in the Stages 
oO . —_|————_ ——__ —_ | 
Output 
Rate gi g2 se) 1 gaa} ga{ 1 2 —l|a-—-lla 
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+ 4 * ry 
N,+3 Not+6 Nat+l2 N4+24 
ft 0 t t, 


Fig. 2—A 1:4 rate-increasing circuit for a 16-bit data block. 


enhancing the rate four times. Figure 2 represents a rate-increasing 
circuit for the 16-bit data block. The incoming polynomial is 


= om aX", (3) 
and 
= X™{X%(a) + aX) + X" (a, + a3X) 
+ Xa, + asX) + +++ Xs + Qy5X)} 
Uy = X™*%*(X% (aq. + aX + aX’ + aX’) 
+ Xa, + a3X + a,X? + a,X*) + --- 
+ Xr. + aisX + aX’ + a,;X°)} 
Ug = XNANANO V8G. + aX + +++ aX’) 
+ X (ag + aX + ayyX” +++ ay5X')} 
Ug = XVANAANOAN GL WAG. tL GX + aad” +++ Gy5X)}. (4) 


The X*° term indicates the delay for the first bit which passes through 
the extra 


Gd 41g ay ain ay, 


periods corresponding to the four stages of the rate-change circuit. 
(See also Appendix B.) Table IT shows a circuit characteristic for an 
n-bit data stream where the rate change is 1:m. The polynomial calcula- 
tions are presented in Appendix B. When it is desired to change the rate 
by a fraction (m,:m.2), a combination of rate-increasing and rate- 
decreasing circuits may be used in series. It is then necessary to choose 
the main clock frequency (f,) which corresponds to the lowest multiple 
of m:, m,, and the incoming data frequency (f;,). The first section of 
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the combined circuit will enhance the frequency from f;, to f, , and the 
second section will reduce from f, to the desired frequency of f;,™m2/m . 
For instance if the input rate is 2400 baud and it is desired to obtain 
3600 baud, then a clock rate of 7200 cycles per second will be necessary. 
The entire circuit for this example is shown in Fig. 3 for a data block 
32 bits long. 


2.3 Fractional Ratio Rate-Changing Circuits 


It is sometimes possible to combine the two independent circuits into 
one circuit and accomplish fractional ratio rate changing. It becomes 
necessary however to perform a special function of delaying selected 
bit positions by a fractional cycle. An example of this circuit is presented 
for a 2:3 rate-change circuit. 


Example 3. A block of data 8 bits long which arrives uniformly during 
12 clock cycles is to be condensed to a data block of 8 clock cycles. 
The incoming data polynomial uy in Fig. 3a is 


i i Ae aX” a aX ae ae 
+ agX° + a,X?”, (5) 
The polynomial u, after delaying alternate bit positions by half a cycle is 
uy, = X*"{ (ay + a,X) + X*(a, + a;X) 
+ X°(as + asX) + X°(a, + a,X)}. (6) 


After the first stage of the circuit in which the gate g, diverts the two 
data positions for three clock cycles into the lower half of the stage and 
the two data positions for the next three cycles into the top half of the 
stage, the polynomial u. may be written as 


Uz = X?{X™*a, + aX) + X™ (Ka, + a3X)) 
+ XX (a, + asX)) + X"(X(ag + a;X))} 
= XM 6g. + aX + aX + azX°) 
+ X%(a_ + asX + aX? + a,X")}, (7) 


TABLE II—CHARACTERISTICS OF A GENERAL RATE-INCREASING CIRCUIT 





Number of Bit Positions Difference in Periods 
Ratio of Input Diverted by Gates in the Stages 
Rate to ————_—|} $j] —__— ___——— 
Output Rate 9: | g2 | gs | cee gt 1;2 {—] a-—-1! a 




















1:m 2° | 21) 22 | --- 2a] k | 2k | — | 2¢°%k 221k 
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and similarly 
Ug = X™*87(KN+7%(G, + aX + aX” + a,X°) 
+ X"**%(a, + asX + aX”? + a,X°*)} 
m= MANA GX + aX” + agX® 
+ a,X* + a;X° + a,X° + a,X"}. (8) 


The magnetic domain technology is especially suitable for creating an 
effect obtained by delaying the alternate bit position by half a clock 
cycle. The peculiarities of the T-bar circuits’ may be exploited in 
achieving this effect. Figure 3b shows two generators, G, and G,, 
excited by the same generating current in the coil. The ‘“‘seed bubbles’ 
(see Ref. 4) of both the generators rotate with the main driving field in 
synchronism. When the coil is excited every 3/2 cycles, G, and G, 
generate bubbles corresponding to the alternate binary bits. The 
distance between the bubbles in section A or section B of the circuit 
will be three T-bar periods. At the junction of A and B, the effect of the 
split Tee will add a quarter period (i.e., X’”*) for the bubbles in section B 
and subtract a quarter period (i.e., X~’”*) for bubbles in section A, and 
the net effect will be to delay the alternate bubble position by half a 
cycle (X’””). The polynomial in section C of the circuit will already 
have a form identical to wu; in (6). 


20 2! 2e 23 24 
t t t ee 
2400baud Ny No N3 Na Ns insuid 





36oobaud Nit! N3+2 N4+4 N5+8 N{+16 


20 2} 22 23 24 


Fig. 83--An example of combination rate-reducing and rate-increasing circuits 
for a 2:3 rate-change circuit for a data block 32 bits long. 
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STAGE | STAGE 2 
SPECIAL ARRANGEMENT 3 CYCLES mieverne 
pee ae Ce ere elt f2 POSITIONS 14 POSITIONS 


POSITIONS BY HALF CYCLE 


\ Ny Ne 


Ny +1 No+2 
3 CYCLES 6 CYCLES 
2 POSITIONS 4 POSITIONS 

Ni 

LEFT eN 

as ‘ “th 

--~-<{ (Gi }--—* GENERAL L”” STAGE 

eer 

Ni+2t-! 


\5 x 2l-' cycLes, 2¢ POSITIONS DIVERTED 


Fig. 8a—A circuit arrangement for a 2:3 ratio rate change for an 8-bit data stream. 
It is possible to add more stages for longer data streams. 


ee cree ee ce ee ee wee ee ee ee 





LAY ET aD 
DU LTLTL FFP) 
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ea ACTING AS X"4 FoR B 








ONE ROTATION 
EVERY X 
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LLL LLL LA 
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SY 


Fig. 3b>—Two generators, G; and G2, excited by the same coil, generate alternate 
bit positions; i.e., G: generates ao , d2 , aa—and G2 generates a , a3 , as—. 
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III. COLLATING AND DISTRIBUTING CIRCUITS 


To avoid undue complications in the polynomial algebra, one example 
of collating circuit and one example of distributing circuit are presented 
in this section. It is possible to generalize the polynomial calculations 
for any general incoming polynomial wp . 


Ezample 4. Four 4-bit inputs from channels 1, 2, 3, and 4 are to be 
collated onto one channel (see Fig. 4). Let the incoming polynomials be 


Uy = AX° + a,X + aX’ + aX, 
Us = A,X* + asX° + agX° + aX’, 
U3; = AgX° 2 QyX° - AyoX”” =e aX", 


and 
U = Q2X” + disx + Qy4X* ++ dig *, (9) 





SECTION 1 

SECTION 2 

SECTION 3 

Ug 

SECTION 4 

Ao a; Ag ag ao 
eS a ee 2 rer 

TIME 
Us 
Ag Ag A991, 
O12 413, Ay4 A15 
Ug 
% Ao A812 ay a13 42 814 a3 ats ao 

| A a4 ae As ag a Ajo a7 Alp Lis 
ae ee TIME 


Fig. 4—A typical 4—-input, 4-bit collating circuit. 
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then each of the sections 1, 2, 3, and 4 must accomplish a rate reduction 
of 1:4. Let 


N; and N2 be the number of periods in top sections of 1, 
Ni and N3 be the number of periods in top sections of 2, 
N? and N3 be the number of periods in top sections of 3, and 
Ni and Né be the number of periods in top sections of 4. 


According to Section II, the lower halves of sections 1, 2, 3, and 4 are 
(Ni + 3), (V2 + 6), and so on. The polynomial u; between the two 
stages of the first section is 


us = X™"*{(ay + a,X*) + X*(a, + aX”’)}, (10) 
and so.on. The polynomials wy , Uio , Wir , ANd Uy2 May be calculated as 
Ug = X™ fg, + a,X*+ a,X*+ aX} 
Uro = XVM 4G + asX* + agX*® + aX} 
iS a yk Hak Fax | (11) 
Uyg = XN Nf A aX) + AyyX” + ay5X}. 
The periods in sections 1, 2, 3, and 4 are chosen such that 
Ni+Ni=Nt+Ni+9 
No + Nz = Ni+N2 + 10 (12) 
N3 + N2 = Ni+ N24 11. 
If the technology for implementation permits N{ + Ni = 0, then 
Uy = X°(ay + a,X* + a,X* + a,X") 
Uro = Xa, + asX* + agX° + aX") 
Uy, = X (ag + agX* + aX + 9,,X"”) 
the = X (ayy + aygX* + ay,X*® + 5X"). 
The output polynomial u,3; can be written as: 
Urs = X° {ao + aX + agX? + a2X?® 
+ a,X* + asX° + aX? + a43X" 
+ ayX® + aX? + aX? + ayx” 
+ ag,X" + a,X™ + a,,X"* + aX}; (14) 


and X° corresponds to the minimum delay in the circuit. It can be seen 


(13) 


RATE-CHANGE CIRCUITS 3187 


that Nt + Nj} need not equal zero. However, equations (12) must be 
satisfied. Under these conditions a fixed delay between the input and 
output results. Further, it is to be noted that equations (12) must be 
chosen with care to ascertain that the output polynomial. u,; has each of 
its term a,X° satisfying the equation 


(minimum delay) + c = b (15) 


where b is the power of X in any one of the terms a,X° in the input 
polynomial. Physically relation (15) implies that no term in the output 
appears before it has been received at the input terminal. In equations 
(12), as b varies between 0 and 15, the corresponding values of c in 
equation (14) satisfy relation (15). 

Other generalizations of this circuit can be readily obtained by 
changing the number of stages (depending upon the number of bits 
in U,, U2 -*: etc.) and the number of inputs to be collated. Reversing 
the direction of propagation yields a reversal of input and output 
polynomials leading to the distributing circuit. 


Example 5. A typical distributing circuit is shown in Fig. 5. The 
input polynomial is 


Uo = Pe + aX + GX” eee ivex (16) 
The binding conditions on this circuit imply that 
Ni+N,+N3 = N, +N. + Ny + 7.* (17) 


The polynomial calculations lead to the output polynomial 
Uz = X {ay + agX + aX? + agX® + agX* + ayy X*? + aya. X° + aX? 
+ a,X* + a,X° + asX"° + aX" + aX” + a, X” 
+ ay3,X™* + a,5X"°}. (18) 
It can be seen that if the direction of propagation is reversed, the 
circuit performs a collating function. 


IV. REVERSING CIRCUITS 


These circuits are capable of changing the order of data bits in an 
incoming polynomial. If the input polynomial is u = 50% a,X; then 
the output of the circuit can be written as u, = X* 5°? a,_:X'’. To 


* Alternatively the lower half of the entire circuit may be designed to have seven 


additional periods if 
Ni+ NA+ N3! =Ni+No24+ Nz. 
Under these conditions the gating at g2,g3 , and gs is somewhat simplified. 
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Fig. 5—An example of a distributing circuit. 


satisfy the physical constraint that no term of the output polynomial 
appears at the output before it is received at the input 


d= n. (19) 


A specific example of reversing an 8-bit data block is presented in 
this section. The generality of this type of circuit is proved in Appendix C. 


Example 6. It is desired to reverse an 8-bit data block represented as 
Up = AX° + aX + aX” + a3X* 
+ a,X* + a3X° 4- agX* + aX". (20) 
After this data block is processed by the first stage in the circuit (Fig. 6) 
Uy = X™{X(a, + aX) + X*(azs + a2X) 
+ X°(as + aX) + X"(az + a¢X)}. 
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1 : 
U3 = X7 3) a7-_X" 
a 





20 2 22 


Fig. 6—An example of a reversing circuit (the general reversing circuits are 
presented in Figs. 8 and 9). 


Further, uv, and us; can be written as 
Ug = K™N Ng. + aX + aX? + aX’) 
+ X"(a, + asX + a,X* + a,X*)} 
Ug = XN NG, + agX + aX? + a,X° 
+ a,X* + a,X* + a,X° + aX’)}. 
If the technology for implementation permits V, + N.-+ N; = 0, then 
us = X'{a, + aX + a;X° + a,X°® 
+ asX* + a,X° + aX° + aX}. (21) 


V. CONCLUSIONS 


The techniques presented in this paper indicate simple combinations 
of gating and propagation functions to yield any general rate-changing, 
collating, distributing, and reversing of data streams. All the necessary 
gates are generally driven by one binary clock diverting 2°, 2, 2”, --- 
binary position (or ‘‘domains’’) from one branch into one or the other 
branches of three branch nodes distributed systematically in the circuit. 
Further, the circuit configurations are optimal within the. constraints 
of the problem to the extent that a minimum number of gates are 
necessary to accomplish any one of the rate changing, distributing, 
collating, or reversing functions. For this general reason the control 
circuitry necessary for implementation of these configurations is reduced 
to a minimum. 

Bubble technology, charge-transfer, or charge-coupled device tech- 
nology each lends itself to the implementation of such circuits. 
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APPENDIX A 


Generalized Calculation for a Rate-Decreasing Circuit 


Let m be the ratio of the incoming data rate to the desired rate and 
n be the number of data bits in the data stream. For initial calculations 
let us choose n to be 2° where a is an integer. The incoming polynomial 
can be represented as 


Up = UX? + aX + aX + aX? +--+ 1X”. 
After the first stage of the general circuit shown in Fig. La, 
Uy = XX fayX° + aX +++ Ay X" 7} 
+X fay Xe Ay X74} 
where n’ = 2° and k = m — 1. Or u, can be written as 
ty = XX (aX? + aX bees dy X" 7) XVM TY 
(dy X°? + Agr g1X foes Oy X”’")}, 
After the second stage of the general circuit 
Ug = K™*N Xi (ayX° + aX +++ Gyr X” ") 
XP Gg XO bee Ayr X” 4} 
+ XO GX fees Oyranr td”) 
a Xe aeeg |X? re 6,1. X"')} 


where n” = 27”, 
After the (a — 1)st stage 


Ua-1 = De eet NX (gk + a,X) 
— x’"(a, + a,X) +++ Dara (ee a An-1X)} 
and after the ath stage 


Ug = tans {aX°+ ax" + xe 


af “ix aes fiack “fe fig 4 Dig iea| 


since m = k + 1 and a = log}. It can be seen that the rate for the u, 
polynomial is (1/m)th the rate of the incoming polynomial u,. When 
the number of data bits n does not equal a number 2%, the general 
polynomial calculations become more complicated, but the difference 
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of elements in the two branches is still the same, and the gates can still 
be driven by binary counters. When the circuit is used on a repetitive 
basis, the gating sequence of the various gates g:, g2 --- g, must be 
altered to adjust for the value of n. One such sequence of a gate operation 
for a 26-bit, 2:1 rate-decreasing circuit (shown in Fig. 7) is presented 
in Table ITI. 


FIRST 4,4,4, NEEM UNIFORM 
a 16 8,8 2 rae — 1 oo 
N,+16 No+8 N3+4 N4at2 Ns +1 
LAST 10 8, a 4,4,4 ao. 1 1 a 
UNIFORM UNIFORM 


Fig. 7—A typical 2:1 rate-decreasing circuit for a 26-bit data block. 


APPENDIX B 


Generalized Calculation for a Rate-Increasing Circuat 


Let m be the rate increase desired and n be the number of bits in the 
data stream. For initial calculation let us choose n as a 2° where a is 
an integer. The incoming data polynomial may be written as 


Uy = AX + aX” + aX” + aX --- aX ™. 
After the first stage of the general rate-increasing circuit (Fig. 8), 
= X™{X"™ "(ay + aX) +X" 
(az + a3X) +++ a?" (a,_g + ay_1X)}. 


Taste JII—Sequence or GATE OPERATION FoR A 26-Bit Data 
STREAM 





Difference in 





Number of Bit Positions Diverted _ Periods in 
by Gates Different Stages 
gi g2 gs 94 Js 1/2)3 [415 
SAA a heey dee so) GV dy oss 


16, 10; | ;8, 8, 8,2; | 4, 4,2; 72.22 427.41; |e] s}4]}2}1 
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Fig. 8—A general reversing circuit for an n-bit data block. 


After the second stage 
Ug = XNtN8 X38" 3g) +t aX + aX’ + a,X*) +: 
EXO gg A dyaX + dpa X” + 1X") 
= xe Peromvg, ta X + aX? + aX!) + 
+ XOVeLaimvg + An-sX + An2X” + A,-1X°)}. 
After the ath stage, : 

Ue = xh Latin ng, + a,X + aX”? + +++ aX" )}. 
When it is necessary to cover the case in which n is not 2* but any given 
number, it is essential to choose an integer number, a, such that 

2° = n. 
The functioning of the rate-increasing circuit can be proved in the 
following way: 


The power of X associated with any term a; in the converted poly- 
nomial u, can be written as the delay d’ which the circuit introduces and 


j=a-1 


d’= YY (m — 1) +4. 


(See the expression for u, in this appendix and assume >—3 N; = 0). 
The first term a, of the incoming polynomial wo must be delayed by 
((n ~— 1)m — (n — 1)) clock cycles to ascertain that ao is (n — 1) clock 
cycles ahead of a,_, , the last term in wu, . After the conversion to u, 
the zth term will have to be z clock cycles behind a, in u,. This leads 
to the total delay d’’ for ith term as 


d’’ = (n — 1I)m — (n — 1)) +7. clock cycles 
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oe Dime) ase 
2° = Dn = 1)-P4 


j=a-l 


d 2m — 1) +3. 


This corresponds to the delay d’ in the polynomial u, . 


APPENDIX C 


Generalized Calculation for a Reversing Circutt 


Consider a block of binary data with n-bit positions and a reversing 
circuit shown in Fig. 8. To simplify the nature of equations let us assume 
that n = 2° where a is an integer. This constraint can, however, be 
easily relaxed in practice as shown in Appendix A 


Up = aX° + a,X' + aX? +--+ aX" " 
u, = X™{(apX? + aX) + (a,X* + agX*) + + ay gX™** + aX” 1} 
uy, = X™{X(a, + aX) + X%(ay + aX) + ++) X™ yy + An-2X)}. 
Similarly 
thay = EEX Mas be dye pX 00+ aX") 
+X daer + On2X +++ Ann X” ')} 
where n’ = 2°’; n’ = n/2 and finally 


Un = yee > Ga (ee + a, ,¢ + WX" ered AX" *)}. 


The order of the polynomial up is reversed. The minimum delay that 
is essential in the circuit is (x — 1) cycles to satisfy the physical con- 
straint that a,-, appear at the output only after it is received by the 
circuit. The delay of >)$ N; depends on the nature of the basic vehicle 
for implementation. With magnetic domain circuits a certain minimum 
is deemed necessary. With charge-transfer or charge-coupled devices 
the value can be made zero. 

When a rate change and reversing are both desired, the number of 
periods in the lower section of any stage 2 should become 


(N; + a’ + a’ ’-k) 


where k = m — 1 and mis the rate change desired. The general circuit 
for a combination reversing-rate change circuit is shown in Fig. 9. The 
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WHERE K’= m+1,N!= 28-1; Usat-1 
DIFFERENCE OF PERIODS IN UTH staGe = 2b+ 2b-1(m-1) 
= 2b-l(mtt) = 2t-tK’ 


Fig. 9—A combined reversing and 1:m rate-change circuit. 


polynomial calculations are very similar to those presented in Appendixes 
A and B. 
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The Overflow Distribution for Constant 


Holding Time 


By P. J. BURKE 
(Manuscript received July 20, 1971) 


An infinite trunk group split into a finite first-choice group and an 
overflow group ts studied. The equilibrium distribution, at an arbitrary 
instant, of the number of busy trunks in the overflow is obtained for the 
case of Povsson input and constant holding time. Some numerical com- 
parisons of variances and distributions for exponential and constant 
holding time are given. The variance of the overflow was found to be always 
the greater for constant holding time, and in the case of one trunk in the 
first-choice group this inequality is proven to be true analytically. In some 
cases studied, the variances differ markedly—by as much as 50 percent. 
Implications of these results for the traffic engineering of overflow groups 
with nonexponential holding time are discussed. 


I. INTRODUCTION 


We consider an infinite trunk group which is split into a finite first- 
choice group and an overflow group. Calls that find all trunks busy in 
the first-choice group are placed on the overflow group. It is assumed 
that the input is Poisson and that the system is in equilibrium. Under 
these conditions the distributions of the number of calls in the total 
group and of that in the first-choice group are known and, for a given 
load, are independent of the holding-time distribution. For the case of 
exponential holding time the distribution of the number of calls in 
existence at an arbitrary instant in the overflow group is also known, 
having been found by Kosten in 1937." The latter distribution, in 
particular its second moment, is basic to the method of engineering 
overflow groups, often called the “equivalent random” method, pio- 
neered by Wilkinson and Bretschneider. (See Ref. 2 for a description of 
this method.) If this distribution were also independent of the holding- 
time distribution, then the equivalent random method could be applied 
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with uniform validity regardless of the underlying holding-time dis- 
tribution. 

It is indeed tempting to speculate that the overflow distribution has 
this independence property, since the number of busy trunks in the 
overflow group is the difference between two random variables: the 
number of busy trunks in the total group and that in the first-choice 
group, each of which is independent of the holding-time distribution. 

Unfortunately the simple example of constant holding time studied 
here proves that the overflow distribution does in fact depend on the 
holding-time distribution. This result is reminiscent of that of Tange 
and Wikell,® who found that the blocking probability in a grading 
depends on the form of the holding-time distribution. In the case of the 
grading the differences between the blocking probabilities for constant 
and exponential holding times are too small to be of practical signific- 
ance. In contrast, the differences between the respective overflow 
distributions, as measured by the relative differences between variances, 
can be quite large, running close to 50 percent in some cases. 

Unpublished work by N. P. Archer indicates, nevertheless, that when 
the equivalent random method is based on constant holding time, the 
results of engineering a first-stage overflow group are essentially the 
same as when the method is based on exponential holding time. That is, 
for a given configuration of first-choice groups and Poisson offered 
loads the equivalent random method, consistently applied, results in 
the same size overflow group for a given loss probability whether the 
holding time is taken to be constant or exponential. Thus the computing 
aids—that is, algorithms, tables, and graphs used in the application of 
the equivalent random method—which are based on exponential- 
holding-time theory can be used without change, as a practical matter, 
for the purpose of engineering a single overflow group even when the 
holding time is constant. A word of caution, however, is in order. It 
should be obvious that the substitution of an exponential for a constant 
holding-time distribution must be done throughout the procedure, in 
the calculation both of the overflow variances and of the size of the 
equivalent group. If exponential theory is used for the overflow vari- 
ances and constant theory for dimensioning the equivalent group, there 
will be a bias toward too few trunks in the overflow group and the 
service may be significantly worse than that aimed for. If the opposite 
error is made and the overflow variances are based on a constant hold- 
ing-time assumption while the exponential charts are used for estimating 
the size of the equivalent group, the overflow group will be over- 
engineered. The latter erroneous result will occur also if constant- 
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holding-time overflow variances are estimated by actual measurements 
of overflow traffic rather than, as is done presently, from theoretical 
considerations only, and if consideration is not given to the effect of the 
holding-time distribution on the overflow variance. 

It is fortunate that the constant-holding-time case is tractable, since 
constant holding-time represents an extreme point in the set of holding- 
time distributions when these are ordered according to their coefficients 
of variation. It is reasonable to conjecture, for example, that an overflow 
variance for any holding-time distribution whose coefficient of variation 
is less than unity will differ less from that for exponential holding time 
than does the overflow variance for constant holding time. Thus, if the 
equivalent random method is applicable in the case of constant holding 
time, it is a fortiort so when the coefficient of variation of the holding 
time distribution is between zero and one. Furthermore, the fact that 
a constant-holding-time equivalent random procedure yields results 
which are almost indistinguishable from those for exponential holding 
time lends support to the conjecture that in the application of the 
procedure the form of the holding-time distribution may be ignored even 
when its coefficient of variation is greater than unity. 

Although the chief purpose of this study was to gain information 
about the extent of the dependence of the parameters of overflow 
distributions on holding-time distributions, the original motivation was 
provided otherwise. In fact the present investigation was sparked by 
the observation that the formula for the decomposition of the variance 
of overflow traffic resulting from a superposition of independent Poisson 
input streams offered to the same first-choice group, derived rigorously 
by A. Descloux for exponential holding time in unpublished work, is 
valid for any holding-time distribution. This formula is 


Var (y;) = pj Var (y) + pil — piEQ), (1) 
where 
y; = the number of calls in the overflow belonging to the 7th stream, 
p; = the proportion of the offered load in the 7th stream, 
and 


y= De Yi - 
The question immediately arises whether the observation that (1) is 
independent of the holding-time distribution has any application, since 
it involves Var (y), which was known heretofore only for the exponential 
case. To settle this question, a characterization of the overflow traffic 
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in nonexponential cases is required. As it turns out, in the constant- 
holding-time case not only the variance but also the distribution itself 
may be found exactly. 

Kosten’s formula for the probability of y trunks busy in the overflow 
for an offered load of a erlangs offered to c trunks may be written, 
after some simplification, as 


a” GS (—a)’ 1 
aoe ———— (2) 
y! j= 9! ytj+tt—1\@: 

| )@ 


4 


pe = 


-¥ 





t=0 


where (c)) = 1 and (c); = c(e — 1) --- (© —~74+ 1). 
The corresponding formulas for the constant case are 


e° act” c a” 
W., Feet Get y ht Mek ey, y > 0, 
“i! (3) 


and 


eo" c (2a)" 2c a” 
W.0 oa c t {s h! + ba hi 


h=ct+l1 


[= (") - ae if Jee Sa |}. (4) 


(A brief outline of the derivation of (3) and (4) is given in Section II. 
Algebraic details of the derivations of these and subsequent formulas 
are given in the Appendix.) 

It might be remarked that the state probability formulas in the 
constant case are the simpler for computational purposes, since they 
involve only finite sums of positive terms. 

With respect to complexity, the formulas for the moments are 
another story. The means, of course, are equal and are given by 


M, = af,.,(a) = - be zy". (5) 


t= 1! 











Although the variance is the second-order moment of direct interest, 
the second factorial moment is equivalent for our purpose and will be 
given here since it is simpler. For exponential service time this moment 
may be written 


rm =a@ri-9 D&+e41].. (6) 
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(From this expression it is immediate that FS” = c¢ + 1fora =c+1, 
a fact that was found useful in checking computer programs.) For 
constant a time the second factorial moment is 





m= cna) Sa a pe + 14g 4 3) — xe) +e > oe 
1=0 a’ 
min [k,c] k 
h=[k-e] + 


In all the numerical cases studied, it turns out that F, > F{”. In the 
case c = 1, the formulas are simple enough to allow an easy analytic 
proof that this relationship is true uniformly in a. The fact that F,> FS” 
implies that the correlation between the number of busy trunks in the 
overflow group and that in the first-choice group is lower in the constant 
than in the exponential case. This fact is perhaps less surprising to the 
intuition than is the result that the overflow is more variable for constant 
than for exponential holding time. 


II. OUTLINE OF THE DERIVATION 


Owing to the constancy of the holding time, taken here and below 
to be of unit length, the calls present in the overflow group at an arbi- 
trary instant, ty) , are precisely those that overflowed during the preceding 
time interval of length unity. The number of calls present in the first- 
choice group at the instant ¢) — 1 is known to have the truncated 
Poisson distribution (also known as Erlang’s first distribution). We 
condition first on the number of calls present in the first-choice group 
at t) — 1. Next, we condition on the number of arrivals during the 
interval [f) — 1, t)), which, of course, has the ordinary Poisson distribu- 
tion. We now observe that the hang-up or departure instants together 
with the arrival instants, as conditioned, are mutually-independently, 
uniformly distributed on the unit interval. This observation enables us 
to complete the calculation by an application of a ballot theorem. 


Let: 
c = number of trunks in the first-choice group 
a = offered load in erlangs 
m=; = probability of z calls on the first-choice group at an arbitrary 
instant 
p; = probability of j arrivals during a unit of time; p; = e ‘a’/j! 
W,, = probability of x calls on the first-choice group and y calls 
on the overflow group at an arbitrary instant 
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w., = probability of y calls on the overflow group at an arbitrary 
instant 
f(y; 7, 7) = probability of y overflows during a unit interval at whose 


initial instant there are 7 calls on the first-choice group and 
during which j new calls arrive. 


Then we may write 


W220 = Dz eS a{(0; i: x) (8) 
7=0 
and 
Wey = Dery 2 TAYE +Y), Yy> 0; (9) 


and, with a reversal of the order of summation, we have 


wo = Dim 21 ps0; 4, i), (10) 
and 
Wy = Dim De phyiid, y> 0. (11) 


The distribution of the number of calls on the first-choice group at an 
arbitrary instant is independent of the holding-time distribution and 
is given by 


m= @'fin| D carn |, 17=0,---,6, (12) 


while the remaining service times of these 7 calls are independently and 
identically distributed according to the equilibrium excess distribution; 
that is, the remaining service time of each call has the distribution 
function 


F(t) = [ isha (13) 


where H(u) is the service-time distribution function (with unit mean). 
The last two results have been published by several authors; an ele- 
mentary proof has been given by L. Takaes,* whose paper includes a 
bibliography on the problem. 

The implication of (13) in the present case is that the remaining 
service times of the calls initially on the first-choice trunks are inde- 
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pendently uniformly distributed on the unit interval. Furthermore, it is 
well-known that when the number of Poisson arrivals during a fixed 
time interval is given, the individual arrival instants are independently 
uniformly distributed over this interval. Hence any specific sequence 
of the 2 departures and j arrivals during the interval [t, — 1, ¢)) has the 
same probability, namely 1/(‘?’). 

The quantity f(0; z, 7) is the probability that a sequence of arrivals 
and departures has the property that at all times the excess of the 
accumulated arrivals over the accumulated departures is strictly less 
than one more than the initial number of idle trunks, i.e., less than 
c — 2+ 1. Thus the problem of calculating f(0; z, 7) can be recognized 
as a “ballot’’? problem. Successive arrivals and departures are called 
“events,” and we denote by a, and 6, , respectively, the accumulated 
number of arrivals and departures at the rth event. With this notation, 


f(0; 7,7) = Pr{a,< 6B, +e-—i+1, r=1,---,¢@+ 7}. (14) 


The required probability is given as a solution to exercise 3 of Chapter 1 
of Ref. 5. We have, with the usual conventions concerning binomial 


( 


= 0, >a 


IA 


Similarly, 

(G3 14) = Pr tos = Bee a Py Plies Late yt bg} 

2 Pre. <p bee ty ft SA te as yO, (18) 
and thus for7 Sc + y, 


ee eee 
e—-ity e—-t+tyt1 


4) =. = ee (17) 
os (oe) 


We are particularly interested in the marginal overflow state probabili- 
ties, namely the quantities denoted w., and w., . These may be put into 
a form suitable for numerical calculation by substituting the values of 
the expressions 7; , p; , and f(y; 7, 7) into (10) and (11). The final results, 
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obtained after a small amount of manipulation, are shown as (38) and (4). 
(See Appendix.) 


III. MOMENTS 


The results (3) and (4) are considerably simpler than their analogs 
for exponential holding time, shown in (2), which was obtained from 
equation (38) of Ref. 1. It does not follow, however, that the moments 
of the distribution defined by (8) and (4) are simpler in form than those 
for exponential holding time. In fact, although in principle each of the 
moments can be written as a finite sum, it seems to be a rather tedious 
task to obtain closed-form expressions for them. Since only the first 
two moments of overflow distributions are presently of practical interest, 
we shall confine our attention to the two lowest-order moments for the 
case of constant holding time. 

The mean of the distribution {w.,; y = 0, 1, ---} is obviously the 
same as that for exponential holding time. Nevertheless, a direct 
calculation of the mean from the expression for the state probabilities 
is useful as a check on the accuracy of the algebraic manipulation 
indicated previously. That is, as a check, the equation 


re De (18) 


y=1 c! i=0 1! 





should be shown to be an identity, as indeed it is. (See Appendix.) 

The calculation of the second factorial moment was done by direct 
summation. The details are given in the Appendix. 

It was remarked above that a simple proof that F. > F3” can be 
given for c = 1. To do this, we observe that by substituting c = 1 
into (7) we obtain | 


F, {a® — a’ + 2a — 2 + 2°}. (19) 


ety ol 
~I+a 
By substituting c = 1 into (6), we obtain 


3 


eae (20) 


Thus we must prove 


q’ — a + 2a — 2+ 2" L a’ 


l+a at ar (21) 
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After cross-multiplying and simplifying, there results, equivalently, 





(2+a)e°>2—a (22) 
or 
a 
Be 8 
é Sa a (23) 
I+, 


which is true for a positive (take logarithms). 


IV. NUMERICAL RESULTS 


A comparison of the overflow variances from first-choice trunk group 
sizes ranging from 1 to 100 is shown in Table I. These numbers were 
calculated on the basis of (6) and (7). The starred entries, for low offered 
loads and high values of c, were calculated separately because of a 
loss of accuracy in (7) at these values. It should be noted that (8), 
having no subtractions, does not suffer from loss of accuracy from this 
cause and hence was used to obtain the variance by direct numerical 
summation. 

The distributions of tails for several trunk-group sizes and offered 
loads are shown in Figs. 1 and 2. It should be observed that even when 
the variances differ noticeably, as for example for c = a, where in both 
cases shown the variance difference is greater than 20 percent, the tails 
distributions are very close at all probabilities of practical interest. 

Since the (almost indistinguishable) curves for c = 10, a = 20 are 
truncated at a high probability level, it may be well to point out that 
they differ negligibly even at low values. Thus the probability of 
exceeding 25 busy trunks in the overflow is 0.0096 for exponential 
and 0.0109 for constant holding time. A negligible difference in the state 
probability distributions, as seen in the case of c = 10, a = 20, was also 
seen in other cases for which a = 2candc > 10. Fore = 20 and a = 40, 
however, a comparison could not be made because of a complete loss 
of accuracy in (2). Since the relative difference in the variances is even 
less in this case than for c small, it is reasonable to conjecture that (8) 
provides an excellent approximation to (2) in this case. This suggests 
that (3) is useful as an approximation to (2) for a > c, precisely in the 
range where (2) is unsatisfactory for numerical computations. 
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TABLE I—COMPARISON OF VARIANCES OF OVERFLOW FOR CONSTANT 
Houpine TIME AND EXPONENTIAL Houpina TIME 





Mean Exp Var Const Var (VD/VE) 





¢ a Ml VE VD VE/Mi VD/M1 —-1 

1 0.40 0.1143 0.1279 0.1331 1.1190 1.1647 0.0408 
0.60 0.2250 0.2575 0.2704 1.1442 1.2017 0.0503 

0.80 0.3556 0.4120 0.4351 1.1587 1.2236 0.0560 

1.00 0.5000 0.5833 0.6179 1.1667 1.2358 0.0592 

1.20 0.6545 0.7661 0.8127 1.1705 1.2416 0.0607 

1.40 0.8167 0.9568 1.0152 1.1716 1.2431 0.0611 

1.60 0.9846 1.1529 1.2228 1.1709 1.2419 0.0606 

1.80 1.1571 1.3529 1.4334 1.1692 1.2387 0.0595 

2.00 1.3333 1.5556 1.6458 1.1667 1.2343 0.0580 

2 0.80 0.1208 0.1478 0.1595 1.2240 1.3211 0.0794 
1.20 0.2959 0.3778 0.4132 1.2767 1.3966 0.0939 

1.60 0.5278 0.6873 0.7562 1.3021 1.4327 0.1003 

2.00 0.8000 1.0489 1.1555 1.3111 1.4443 0.1016 

2.40 1.1006 1.4425 1.5864 1.3106 1.4413 0.0997 

2.80 1.4218 1.8551 2.0330 1.3048 1.4299 0.0959 

3.20 1.7579 2.2784 2.4858 1.2961 1.4140 0.0910 

3.60 2.1054 2.7075 2.9392 1.2860 1.3960 0.0856 

4.00 2.4615 3.1392 3.3902 1.2753 1.3773 0.0800 

5 2.00 0.0734 0.1040 0.1198 1.4176 1.6319 0.1512 
3.00 0.3302 0.5186 0.6078 1.5707 1.8409 0.1720 

4.00 0.7963 1.3013 1.5284 1.6342 1.9194 0.1745 

5.00 1.4243 2.3332 2.7210 1.6381 1.9104 0.1662 

6.00 2.1624 3.4864 4.0166 1.6123 1.8575 0.1521 

7.00 2.9730 4.6819 5.3164 1.5748 1.7882 0.1355 

8.00 3.8321 5.8806 6.5798 1.5346 1.7170 0.1189 

9.00 4.7242 7.0661 7.7970 1.4957 1.6504 0.1034 
10.00 5.6395 8.2327 8.9712 1.4598 1.5908 0.0897 

10 4.00 0.0212 0.0329 0.0399 1.5485 1.8816 0.2151 
6.00 0.2589 0.4872 0.6044 1.8821 2.3349 0.2406 

8.00 0.9733 1.9857 2.4530 2.0402 = 2.5203. 0.2354 
10.00 2.1458 4.3624 5.2907 2.0330 2.4656 0.2128 
12.00 3.6231 7.0710 8.3540 1.9516 2.3058 0.1815 
14.00 5.2820 9.7876 11.2465 1.8530 2.1292 0.1491 
16.00 7.0490 12.4046 13.8965 1.7598 1.9714 0.1203 
18.00 8.8826 14.9080 16.3499 1.6783 1.8407 0.0967 
20.00 10.7592 17.3130 18.6671 1.6091 1.7350 0.0782 

20 8.00 0.0013 0.0021 0.0027 1.6141 2.1127 0.3089 
12.00 0.1175 0.2584 0.3401 2.1986 2.8935 0.3161 
16.00 1.0306 2.7028 3.5055 2.6226 3.4015 0.2970 
20.00 3.1778 8.2920 10.3872 2.6093 3.2686 0.2527 
24.00 6.1700 14.8143 17.7126 2.4010 2.8708 0.1956 
28.00 9.5977 20.9311 23.9282 2.1808 2.4931 0.1432 
32.00 13.2464 26.4675 29.2087 1.9981 2.2050 0.1036 
36.00 17.0146 31.5567 33.9567 1.8547 1.9957 0.0761 
40.00 20.8522 36.3397 38.4213 1.7427 1.8426 0.0573 
50 20.00 0.0000 0.0000 0.0000* 1.6452 2.2408*  0.3403* 
30.00 0.0066 0.0161 0.0230 2.4215 3.4727 0.4341 
40.00 0.7476 2.7343 3.7810 3.6573 5.0574 0.3828 
50.00 5.2393 19.7746 25.5760 3.7743 4.8817 0.2934 
60.00 12.9671 40.9397 48.5749 3.1572 3.7460 0.1865 
70.00 21.9661 57.8384 64.1752 2.6331 2.9216 0.1096 
80.00 31.4446 71.6558 76.4680 2.2788 2.4318 0.0672 
90.00 41.1339 83.8760 87.5797 2.0391 2.1291 0.0442 
100.00 50.9303 95.2624 98.2036 1.8704 1.9282 0.0309 
100 60.00 0.0000 0.0001 0.0001* 2.4634 3.6885* 0.4973* 
80.00 0.3194 1.4158 2.0681 4.4331 6.4758 0.4608 
100.00 7.5700 38.5958 50.7901 5.0985 6.7094 0.3160 
120.00 23.5523 89.6717 103.9455 3.8073 4.4134 0.1592 
140.00 42.1741 123.6750 133.0991 2.9325 3.1559 0.0762 
160.00 61.5406 149.8224 156.0855 2.4345 2.5363 0.0418 
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Fig. 1—Tails distributions, c = 10. 
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Y= NUMBER OF BUSY TRUNKS IN OVERFLOW 


Fig. 2—Tails distributions, c = 20. 
APPENDIX 
A.1 Formulas for the Overflow State Probabilities 


The probability of zero calls in the overflow group is given by (10), 
which, after substitution of the values of 7; and p,; , becomes 


c z t c 7 
Wo= dD oe 2d e* 10; 1, j). (24) 


7=0 





h=0 h! 
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For brevity, let 
¢ kh |-l 
f(a) = | > | : (25) 


Then, replacing f(0; 7, 7) by its value as given in (15), we may write 





¢ ait 
W.0 3 f(a) Pp» a! j! 


(#41) 


Be cs | | eee | 
- 10 S53 creer = 


Setting ¢ + 7 = h and multiplying the numerator and denominator 
of each term by h!, 


a” h h! 
wo = Ko] O¢ hi ath oo! ih — 1)! 
2c a” Al h! 
uv a hl <4 ere (h—a! (h gras tae (27) 


w.0 = f(a) Sy & iL, (3) = Ge b+ Wl, | 


h=ct+1 i=h-c 
(28) 


a 
c 


Use of (25) now yields (4). 

The probability of y calls in the overflow group is given by (11). 
We first simplify f(y; 2, 7), as given by (17). We have 
f¢ “a )) als Serer eT) | ee ee 

WOU T+I- CHOC +H! 

a! 7! (29) 
are (eer ee Diep yp ly! 

Thus, from (11), 





eS a/tl SS et! a!9! | 
a = > a’ j-A%, 2LeE+I-CH+YDI!C+yYy! 
tao A} 


eee ek ya DI Me-P yb iy (30) 


& 
I 


c+y ar 1 
f@ & da lex3 SC re ra! 


1 
sae 


Setting h = 271+ j—-— (c+ y), reversing the order of summation, and 


OVERFLOW DISTRIBUTION 3207 


factoring yields 


gt 


vw, =10 D Daegyaepietyti-o, (32) 


h=0 i=h 


ct+y 


a a" 
w» =1% CTF DI = Me+y+ 1—h), — (33) 
and (3) follows. 


A.2 Proof That (18) Is an Identity 
Let the mean of the overflow distribution be denoted by M. 
From (3), 


fe) 


MI) = Dee LHe it Metyt1-) GH 


y 


ty 


M/f(@) = 2s ee oN oe 


M/K(@) = bye 1) bet 


etyt1 


ya 
-~XG gai e - pT) Darren (36) 


After changing the indices in the subtracted sums by replacing 
jby7+ 1 and y by y — 1, we obtain 


+y 


M/f(a) = > oo — Pad) hero — 


a’ : ao 
— —(¢ — — 1) ————- 37 
70 1 DL erp! (37) 
Simplification yields 


MHN@) = os > it Fe Isr 


Putting h = y + cand adding and subtracting terms, we have 
c a’ ) a" c “| 
M/f(a) = RilEe-9h- Bo aA 
i) c h 
+ EFe- fee - |. (39) 
pear Le rs OF 
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We notice that the finite series vanish and the infinite series combine 
to give 











c ait! c-1 | 00 a* att! : 
fe Pp eee a 


Replacing f(a) by its explicit expression, we have 


wee (= ay" (41) 


= 8! 





A.3 Derivation of the Formula for the Second Factorial Moment of the 
Constant-H olding-Time Overflow Distribution 


From (3), 


ety 


F,/f(a) = > yy — 1) Gey i 


Fetyt1-Ie-F+) 42) 


ety 


Poa) = UEe-i+) Yww-Ye+y+V eq yaw 


ct+y 


— ye (-j+) > uy - i eaeereriary (a (43) 


Again replacing 7 by 7 + 1 and y by y — 1 in the subtracted sums 
in (43), we obtain 


ct+y 


F,/f(a) = SFe-it) Dv - DeareesT 


ety 


= Bae Su-Do-Iey_CiCD 


exif 
(c+ y)! 
+20 5¢-)Lo-Ve yw) 


Letting h = c + y and adding and subtracting terms allows us to write 


Pi = S2{S @-on-0- 9% 
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a” <a 
74h —oh—e-Ip+2 de — 5) 
ie) h c h 
4D a@-e-YS- Vdc) (46) 
h=0 h! h=0 hi 


F./f(a) = ys [h(h — 1) — 2he +c + ms rs 


hoc 7 


+ > ale — oe + IF > 1 
— 2, Pak — Dao + 1 2 


- dF “Do Se- e+ 1 +h) 


+2 %@-9 ie tl— 4. (47) 


After substituting & for h + j in the finite sums and substituting e* for 
its series expansion, we obtain 


F./f(a) = e{[2a( +1)-—a@—ce+))] 


at [2a° — 2a(e + 1)] a“ 


7=0 J 


2c-1 a® min [k,c] 


k 
+ Da —> ( Je+1 —h)\(c — 2k + 3h). — (48) 
k! h=[(k-ec]+ h 
It turns out that the second sum vanishes for k = 2c — 1, and thus 
after deleting the vanishing term and writing the coefficient of e* as a 
polynomial in a, we have 


e+2 e+l1 


a 


aj) 
F,/fla) = e{ ek (c — 1)! 


— DF WG-1) - ze +) +ee+ vi} 








2c-2 a’ min [k,c] 


ee ds (e+ 1 mle — 2k + 3h). (49) 


h=[k- 


Finally, after substituting the value of f(a), given by (25), transposing 
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the limits of summation in the first sum ‘on the right-hand side, and 
some minor simplification, we obtain (7). 
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Statistical Behavior of a Fading Signal 


By S. H. LIN 
(Manuscript received July 22, 1971) 


A general analysis of the statistical behavior of the envelope of a fading 
signal V(t)e’*” ts presented in this paper. The statistics include the 
probability P(V < L) that the amplitude V(t) will fade below a specified 
signal level L; the expected number N(L) of fades of V(t) below L per unit 
time; and the average duration t(L) of fades below L. The model for the 
fading signal zs a constant vector plus a random interfering vector which 
represents the resultant of all the received extraneous signals and noise. 
The theoretical results agree with three empirically observed power relation- 
ships obtained in deep fades of nondiversity signals: PV S L) « L’, 
N(L) « Land t(L) « L. The theoretical results are applicable to a wide 
class of fading problems. The analysis includes the previous works of Rice, 
Nakagamt, Norton, Vogler, Mansfield, and Short as special cases. 


I. INTRODUCTION 


A general analysis of the statistical behavior of the envelope V (t)e’*“” 
of a fading signal is presented in this paper. Our principal interests are 
the probability, P(V < L), that the amplitude V(¢) will fade below a 
specified signal level L* ; the expected number, N(L), of fades per unit 
time below the specified level L; and the average duration, f(L), of 
fades below L. These statistics are all functions of the signal level L. 

The theory presented herein has been developed to complement and 
extend the empirical results developed by my colleagues at Bell Tele- 
phone Laboratories’ ° from their extensive experimental experience. 
Published data of other workers have also been considered. 

The previous theoretical works on the statistics of a fading signal 
often assume a complex Gaussian model for the fading signal. The 
theoretical support for this assumption is that, by the central limit 
theorem, the real and the imaginary parts of the sum of a large number 
of independent interfering signals will be approximately Gaussian. 


* More precisely speaking, in a long time period containing a large number of 
fades, the distribution P(V S L) represents the expected fraction of this long time 
period that the signal amplitude V will fade below L. 
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For tropospheric radio links, this model seems to be satisfactory. 
However, for line-of-sight radio links, the results of a short pulse 
experiment’ and the angle-of-arrival measurements’’’’’® indicate that 
the number of interfering signals is usually fairly small. Ray tracing 
theory also indicates that for typical line-of-sight radio links the number 
of paths contributing to multipath propagation is unlikely to be large. 
Furthermore, the theoretical results of the complex Gaussian model 
do not agree well with the experimental data on the statistics of fading 
signals of line-of-sight radio links, especially for certain overwater paths 
with severe fading. 

In this paper, we do not impose the restrictive assumptions of the 
complex Gaussian model. Rather, we simply model the fading signal 
Ve'* asa constant vector plus an interfering random vector; i.e., 


Ve = 14+ Re =1+a+ 38, (1) 


where FR, 6, a, and @ are the amplitude, phase, real part, and imaginary 
part respectively of the interfering vector. The interfering vector is 
described by the joint probability density function f(a, 8) and represents 
the resultant of all the received extraneous signals, echoes, rays, and 
noise. The analysis applies for R and @ either dependent or independent; 
6 uniformly or nonuniformly distributed; a and 6 either Gaussian or 
nonGaussian. Thus, the results of this analysis may be applied to a 
wide class of fading problems. 

This paper treats the problem in three parts: The first is concerned 
with the amplitude distribution of V. The second considers the number 
of fades N(L) and the average fade duration f(L). The final section 
investigates several special topics including m-distributions, chi- 
distributions, the Rayleigh distribution, Rice distribution, log-normal 
distribution and the sum of n unit random vectors. 

Appendix A is a list of symbols and their definitions. 


Il. SUMMARY OF RESULTS 


(t) In spite of great variations in fading environment and test 
conditions, the experimental data’°"*~**""°"" on P(V S L), 
N(L), and i(Z) of most nondiversity* fading signals obey the 
following three prevailing power laws of deep fades: 


* The diversity signal is the output signal of a diversity combining system with 
two or more input signals. A “nondiversity fading signal” is a fading signal that is 
not a diversity signal. 
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PVsDe«eL (2) 
N(L) « Lf forsmall L. (3) 
i(L) « “| (4) 


The theoretical analysis shows that if the probability density 
function f(a, 8) of the resultant interfering vector, Re’? = a + jp, 
is a smooth function which is neither singular nor zero at the 
deep fade point (a2 = —1, 8 = 0), then the statistics P(V < L), 
N(L), and ¢(L) of deep fades follow the three prevailing power 
laws (2), (8), and (4). The easily satisfied condition, © > 
f(—1, 0) > 0, is sufficient to obtain these functional relationships. 

(77) The set of power laws (2), (8), and (4) apply for R and 6 either 
independent or dependent, @ either uniformly or nonuniformly 
distributed, (or, a and 8 either Gaussian or nonGaussian, either 
independent or dependent) as long as f(a, 8) is smooth. 

(iz) If f(a, 8) is singular at (a = —1, 6 = 0), then the theory predicts 
that for small Z 


P(V SL) « L™ | | (5) 
N(L) « D**e 1>y>h. (6) 
i(L) « L (7) 


The exceptional behavior (5) consistent with » = 1/2 has been 
observed experimentally on certain overwater radio links with 
severe fading. In this case, the resultant interfering vector 
contains the strong water-reflected ray as a dominant compo- 
nent. Therefore, the probability density function f(a, 8) has a 
singularity at the position of the dominant component vector.* 


(iv) If f(a, B) has a zero at (a = — 1, 8 = O) or is negligibly small 
at (a = —1, 6 = 0), the theory predicts that for small L: 

P(V SL) « L* | (8) 

NIL) « 1 p> (9) 

UL) « L (10) 


The composite fading signals of diversity combining systems 
obey the set of power laws (8), (9), and (10). For overland radio 


* The complex Gaussian model, which assumes that f(a, 8) is a two-dimensional 
normal density function, is unable to explain the exceptional behavior (5) and (6). 
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links, the parameter yu is equal to the order of diversity. In our 
experiments, the nondiversity fading signal of a relatively short* 
radio link with path length 15.87 miles has also shown the 
exceptional behavior described by (8), (9), and (10). 

(v) The theoretical results (4), (7), and (10) indicate that the power 
law, i(L) « L, for the average fade duration is more universal 
than those of P(V < L) and N(ZL). This prediction agrees with 
available experimental data. 

(vt) In general, the relationship between f(a, 8) and the amplitude 
distribution P(V S JL) is not unique. As an example, this non- 
uniqueness shows that specifying the Rayleigh distribution for 
the amplitude of a fading signal does not necessarily imply that 
a and 6 are Gaussian, nor does it necessarily imply a large 
number of interfering signals. 


* At 4-GHz operating frequency, the average path length of line-of-sight radio 
links is about 27 miles. 
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Part 1. Amplitude Distribution 
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I. INTRODUCTION AND SUMMARY 


In the study of fading signals due to multipath interference, the 
experimental data on the cumulative amplitude distribution, P(V S L), 
of deep fades are often plotted on a graph paper where the fade depth 
is expressed in dB and the probability of fade is expressed on a log scale 
as shown in Fig. 1. The consensus based on large amounts of experi- 
mental data’~°'**~**'°"”” is that the cumulative amplitude distribution 
of most nondiversity fading signals in the deep-fade region can well be 
represented by a straight line with a prevailing inverse slope of 10 dB 
per decade of probability. 

The equation which describes this typical distribution on Fig. 1 is 


PV <L)=e¢L’, for L, =>L=0 (11) 


where V is the envelope voltage of the random fading signal normalized 
to its nonfaded signal level, L is any specified signal level, ¢ is a param- 
eter depending on fading environment, and J, is the upper bound of 
signal level below which the straight-line representation of P(V S ZL) on 
Fig. 1 is valid. 
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Fig. 1—Cumulative amplitude distributions of fading signals. 


The empirical result (11) means that the amplitude distributions of 
most nondiversity fading signals obey the following square law of deep 
fades 


PVSI)«L’, L,~zZLz20 (12) 


in spite of the great variations of fading environment and test condition. 
However, there are some exceptional cases. The most profound excep- 
tion occurs on certain overwater radio links with severe fading. In 
these instances, the probability of fades, P(V < L), decreases very 
slowly as the signal level L decreases and is characterized by an inverse 
slope of 20 dB per decade of probability as shown in Fig. 1, implying 
a power law 

PV sD «QL (13) 


in the deep-fade region. 

Another kind of exception occurs on certain radio links with relatively 
little multipath fading. The probability of fades, P(V S L), decreases 
very rapidly with Z and is characterized by an inverse slope of 5 dB 
per decade of probability as shown in Fig. 1. This kind of distribution 
follows the power law 


PV sL«L' (14) 
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in the deep-fade region. For example, at 4-GHz operating frequency, 
the average path length of line-of-sight radio links is about 27 miles. 
The behavior (14) has been observed to occur on a relatively short path 
with path length of 15.87 miles. 

The theoretical amplitude distributions previously derived were based 
on a complex Gaussian model and predict a square-law dependence, 
P(V S L) « L’, in the deep-fade region. For example, the Nakagami 
distributions,”® which include Hoyt distribution,”’** Rice distribution,” 
and Rayleigh distribution as special cases, all are square law in the 
deep-fade region. The explicit expression of Nakagami distribution can 
be found in Equation (4.6-28) of Reference 18. Figure 1 also includes 
Rayleigh distribution (dashed line) for comparison with the experimental 
data. 

The small number of interfering signals in line-of-sight radio links 
suggests that the assumption of a complex Gaussian model may be 
unjustified. One of the main objectives of this paper is to determine the 
weakest set of assumptions under which the square law (12) is obtained, 
and the condition for which the exceptional case such as (13) or (14) 
will occur. 


In summary: 


(7) The theoretical model for the fading signal, Vef, is a constant 
unit vector plus a resultant interfering vector 


Vel@ =1+ Re”? =1Llt+at jp. 


The resultant interfering vector Re’? = a + j8, with joint 
probability density function f(a, 8), represents the sum of all 
the received extraneous signals, echoes, rays, and noise. 


(iz) An infinite fade (i.e., V = 0) occurs whenever a = —1land 8 = 0. 
At this point R = 1 and @ = zw. Therefore, the behavior of 
f(a, 8) near the infinite fade point (a = —1, B = 0) is closely 


related to the power law of amplitude distribution P(V S L) 
in the deep-fade region. . 

(zit) For most radio links, the interfering signals and noise may be 
considered random, so that the joint probability density function 
f(a, B) of the resultant interfering vector is a smooth function 
near the infinite fade point (a = —1, B = 0). The analysis 
shows that if f(a, 8) is a smooth function which is neither 
singular nor zero at (a = —1, 6 = 0), then P(V S L) « L’ for 
small LZ. The simple condition » > f(—1, 0) > O is easily 
satisfied by most radio links. The validity of this square law does 
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not require that f(a, 8) be a normal density function. Therefore, 
the number of interfering signals does not have to be large. 

(iv) The analysis shows that if f(a, 8) is not smooth, but is singular 
at (a = —1, 8 = 0), then for small L, 


PV SL «L™ 1>pn2h. 


A physical example for this case is the overwater radio link 
where the water-reflected ray is almost as stable as the direct 
ray. The resultant interfering vector in this case contains the 
water-reflected ray which is “not very random.” In other words, 
the joint probability density function f(a, 6) has a high peak 
at the position of the dominant, stable component, and. may be 
considered singular at that point. 

(v) If f(a, B) is zero or is negligibly small at the infinite fade point 
(a = —1, 8 = 0), then the analysis shows that for small D 


PV SL «<L™, p>. 


A physical example for this case is the short radio link where the 
phase differences among the multipath propagations are all 
small. Then the value of f(—1, 0) is negligibly small because 
the phase, 6, of the resultant interfering vector is generally 
small. Another example for this case is the composite signal of 
the output of a diversity combining system where the artificial 
active combining device serves to create a zero of f(a, B) at 
(a2 = —1,B = 0). . . 


II, FADING SIGNAL MODEL 


The received fading signal is modeled as a constant vector plus an 
interfering random vector as shown in Fig. 2. The latter represents the 
resultant of all the received extraneous signals, echoes, rays, and noise. 
The received fading signal normalized to the magnitude of the constant 
vector can be written as | 


Ve? =j+ Re, (15) 


where R(t) and 6(¢) are the normalized magnitude and the phase of the 
interfering random vector respectively; V(t) and ¢(é) are the normalized 
magnitude and the phase of the received fading signal respectively. 

Let x(t) and y(t) be the real part and the imaginary part of the com- 
plex fading signal Ve’*, i.e., 


Ver = x) + jy; (16) 
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and let a(t) and B(t) be the real part and the imaginary part of the 
complex interfering random vector Re’’, i.e., 


Rite’? = a(t) + JB). (17) 


V,o, R, 6, x, y, a and @ are all real random variables. The normalized 
output power is 


=1+ 2Rcs0+ RR’ = (1+ a)’ + 8B’. (18) 


The relative phase, 6(¢), between the interfering vector and the constant 
vector can be taken to have values from 0 to 27 because @ and (@ + 2nz) 
for any integer n are indistinguishable to the received signal at an 
operating frequency. 

A geometrical interpretation of equations (15) to (18) shows that 
deep fades (i.e., small V) occur when R and 6 are near the infinite fade 
point, (1, 7), in the (R, @) plane or equivalently when a and 6 are near 
the infinite fade point, (—1, 0), in the (a, 8) plane. 

For line-of-sight radio links, notice that as far as the received signal 
V(ée’* is concerned, scintillation, atmospheric divergence, and earth 
bulge effects may also be replaced by a mathematically equivalent 
interfering signal R(t)e"’“ which when combined with the constant 
vector gives the received fluctuating signal. 

_ Therefore, the model described by equations (15) to (18) alae: 
many possible fading mechanisms which may occur individually or 
simultaneously on a fading environment. 

Going one step further, this model also includes the situations where 
no physical constant vector* exists because the mathematical decom- 
position (15) is applicable to any arbitrary fading signal. For such 
situations, the constant vector may represent the average signal level. 
The fluctuation of the signal is considered to be caused by an equivalent 
resultant interfering vector Re’. 

For radio links subjected to multipath interference, the number of 
incoming component waves is usually more than two. Thus the resultant 
interfering vector Re’’ consists of more than one echo and should not be 
interpreted as a simple physical echo. The main distinction is that the 
magnitude and the time delay of a physical echo are not functions of 
operating frequency whereas the magnitude R and the equivalent time 
delay, ¢, = 0/w, of a resultant interfering vector consisting of more than 
one echo, are functions of operating frequency (i.e., are dispersive). 


* For beyond-the-horizon radio links, there is no direct radio path between the 
transmitter and the receiver. 
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Fig. 2—Fading signal model. 


Ill. GENERAL FORMULATION OF AMPLITUDE DISTRIBUTION 


Equation (18) shows that the probability that the random signal V 
be faded below a specified signal level L is equal to the probability that 
a and @ fall‘within the circular region 


(+e’?+e@esLD (19) 


in the (a, 8) plane as shown in Fig. 3. Let f(a, 8) be the joint probability 
function of a and 8. Then P(V S L) is the integral of f(a, 6) over the 
circular region (19); ie., 


B=L a=-1+VL2—p? 
PIV <1) = i _ [ pee Hes 8) da ag. (20) 





The statistical behavior of the interfering vector, Re” = a + j8, is 
sometimes described by the joint probability density function q(R, @) 
of the magnitude R and the phase @ of the interfering vector. A similar 
derivation in terms of R and @ yields 


R=1+L G=r+O4n 
PV <b = i i o(R, 6) dé aR, (21) 
R=1-L O=7—-OL 
where 
syste lee E) 
6, = cos QtH AE . (22) 


Most of the following analysis is in terms of a and 8. An equivalent 
result in terms of R& and @ is given in Appendix B. 





Fig. 3—Fading signal model on (a, 8) plane. 
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IV. PROBABILITY DENSITY FUNCTION OF INTERFERING VECTOR 


Equation (18) indicates that the infinite fade (i.e, V = 0) occurs 
when a = —1, 6 = O. For most overland radio links, the interfering 
signals and fhe noise vary continuously in a random manner, so that 
the joint probability density function f(a, 8) of the resultant vector is 
a smooth function. On the other hand, there are some paths for which 
f(z, 8) may not be smooth, but singular. For example, for overwater 
radio links, the water-reflected ray near the grazing angle is comparable 
in magnitude and stability to the direct ray. Therefore, the resultant 
interfering vector contains a dominant component and the probability 
density function f(a, 8) has a sharp peak (i.e., singular) at this point. 
If the heights of the antennas and the path length are such that the 
average phase of the water-reflected ray is equal to 7, then the singu- 
larity of f(a, 8) occurs at the infinite fade point (2 = —1, B = 0). 

These discussions suggest that in the general analysis of deep fades, 
one should consider not only the case with smooth f(a, 8) but also the 
case where f(a, 8) is singular at (a = —1, 8 = 0). A general probability 
density function f(a, 8), which is useful in our study, is*t 


f(a, 8) = [(l + a)’ + B")""-Ha, 8) (23) 
= Ve"). Ala, 6). (24) 


where H (a, 8) is an arbitrary smooth function. In the range 1 > p 2 1/2, 
the density function f(a, 8) has a singularity of order 2|y — 1] at 


(a = =1,p = 0). 

On the other hand, in the range © > yw > 1, the density function 
has a zero of order 2(u — 1) at (a = —1, 6 = 0). For convenience, we 
shall call the parameter, u, the smoothness index of f(a, 8). 

Since the possible singularity or zero of f(a, B) at (2 = —1, B = 0) 
is taken care of by the factor [(1 + a)” + 67]*"', we shall assume that 


o > H(—1,0) > 0. (25) 


Thus, the density function f(a, 8) given by equation (23) is neither 
singular nor zero at (a = —1, 8 = 0) if, and only if, » = 1 because the 
condition » = 1 implies f(a, 8) = H(a, B), and vice versa. 

In equations (23) and (24), the smoothness index p» can be either an 


* The reason for the use of the factor (u — 1) instead of a simple power index in 
equations (23) and (24) is for the convenience of notation in Section IV of Part 3 
when we investigate the m-distributions. 

{t The cases where the marrage of f(a, 6) occurs at positions other than the 
infinite fade point (a = —1, 8 = 0) will not be analyzed in this paper. A brief 
discussion is included in Sei VII of this part. 
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integer or a noninteger. The only restriction on yp is that 
p= 3%. (26) 


The reason for this constraint on p» is that in Section V of Part 2 we 
find that if u < 1/2, then the expected number of fades V(L) approaches 
infinity as the fade depth ZL approaches zero. This seems to be non- 
physical. Therefore, we require that » 2 1/2.* 

Since H(a, B) is an arbitrary smooth function and uw can range from 
1/2 to «, then the probability density function f(a, 8), as given by 
equation (23), includes a large variety of fading environments. 


V. POWER SERIES REPRESENTATION OF AMPLITUDE DISTRIBUTION 


We shall assume that H(a, 8) is sufficiently smooth so that the two- 
dimensional Taylor series” expansion of H(a, 8) is applicable in the 
neighborhood of a = —1, and 6 = O. Several situations, where the 
Taylor series expansion of H(a, 8) is not applicable, will be discussed 
in Section IX and Appendix C. 

The Taylor series” expansion of H(a, 8) gives 


n=0 


He, 6) = "Se | 4, ¥ crtt..(-1, OC + a)'8 | 


=0 


A(—1, 0) al H, (—1, 0)(1 i a) Se A, .(-1, 0)8 (27) 


al i (A. ,o.(—1, 0) + a)” + 2H,.(—1, 0)(1 +a)6 + Hy,.(—1, 0)67] 


where 
Hy-r(1, 0) = ser gg He, B)| (29) 
B=0 
eo. n! 
C: ~ rl(n—r)t (80) 


Substituting (27) and (23) into (20) for P(V S LZ) and carrying out the 
integration (Appendix D) yields 


PV s -> deel (31) 


* Notice that the unity total probability requires that the singularity of f(a, 8) 
be integrable (i.e., 1 > 0). The constraint 1 = 4 does not violate this condition. 
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= d,L”* ae > deals?" (32) 
S=1 
where 
pox awH(—1, 0) (33) 
W 
ms = Hy s-2¥.2,(—1, 0) 


(34) 


des+2 = 


SiS + w)2% oo |=) GS — »)! 
S! = SCS — 1)(S — 2). . 8-2-1. 
The corresponding amplitude probability density function p(L) is 


S=0 
pL) = 2 PW sD) = YS +2) dasyal? 7". (85) 


In the deep-fade region where LZ is small, the limiting forms of 
P(V s L) and p(Z) are 


L-0 L 
and 
p(L) —> 2rH(—1, 0)L**", (37) 
L-0 


where » 2 1/2. The power law of deep fades for the three different 
cases (2) wp = 1, (4) 1 > w = 1/2, and (iz) w» > 1 with their physical 
fading environments will be discussed in the following sections (VI, VII, 


and VIII respectively). 


VI. PREVAILING SQUARE LAW OF DEEP FADES 


For the nondiversity fading signals of most radio links, f(a, 8) is 
neither singular nor zero at (a = —1, 8 = 0). Theny = 1 and f(a, 8) = 
H(a, 8). Equations (31) to (37) under this situation become 


PVE = Stes (38) 
=rf(-1or+daLl'+daL°+::-, (89) 

f(a, B) = H(a, 8), (40) 

d, = mC=i, 0) = rH (—1, 0), (41) 

ee = Tv > (Pre re eGo 0) (42) 


GS + 1)! on v=0" vN)(S 5 y)! : 
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are 


fes—2v,2(—1, 0) = DS BS -a7 nAay f(a, B) ’ (43) 
0a op v3 ay 
p(L) a D> (2S + 2) dsgval (44) 
S=0 
PV s L) —af(-1, OL’, (45) 
L-0 
and 
p(L) ——> 2rf(—1, O)L. (46) 
L-0 
Equation (45) means that as long as f(a, 8) is neither singular nor zero 
at (a = —1, B = 0), then the cumulative amplitude distribution in the 
deep-fade region always obeys the square law 
PVSDI)«L’,, Li, 2zLezo. (12) 


Notice that this conclusion does not depend on any specific probability 
density function f(a, 8) for the interfering vector as long as f(a, 8) is 
smooth and © > f(—1, 0) > 0. The conclusion applies for a and 6 
either normal or not, either dependent* or independent, either with 
zero mean or with nonzero means. The magnitude FR and the phase @ of 
the interfering vector can be either dependent* or independent and 
6 can be either uniformly or nonuniformly distributed. Therefore, this 
conclusion covers a wide class of signal fading problems. 

Apparently the simple condition, « > f(—1,0) > 0, is appropriate to 
the nondiversity fading signals of most radio links because the square law 
of deep fades is representative of the experimental data.*~°°?~*4:7°"**t 

Notice that the first terms of equations (39) and (45), wL’, are the 
area of the two-dimensional region on (a, 8) plane bounded by the 
circle, L? = (1 + a)? + 8’, in which V < L, as shown in Fig. 3. 

The coefficient f(—1, 0) in equation (45) has been observed to depend 
upon path length, operating frequency, path profile, and geographical 
factors. From the experimental data of a large number of radio links, 
it is possible to deduce an empirical formula of f(—1, 0) as a function of 
these parameters.77'7°'74:?° 

In equation (12), the upper bound, Z,, , of signal level below which 

* This conclusion does not hold if the correlation coefficient between a and 8 or 
between R and @ is unity because the joint probability density f(a, 8) becomes 
Ot This theoretical result also explains an experimental fact that the observed 


amplitude distributions of atmospheric radio noise are also characterized by the 
square law (12) in the small amplitude region.! 
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the square law applies also depends upon fading environment. Our 
experimental data show that L,, of most line-of-sight microwave radio 
links is above 0.3 (i.e., above —10 dB). 

If f(—1, 0) is negligibly small, so that the first term in the power 
series (39) can be neglected, then the second term, d,L*, dominates, 
with the result that the amplitude distribution follows the power law 
P(V < L) « L*. We have observed this behavior on a short radio link 
in the signal range from —10 dB to —20 dB (ie, 0.3 = L 2 0.1). 
For fade depths deeper than —20 dB (i.e., 0.1 > LZ 2 0) the quadratic 
term of (39) again dominates; and the transition region between 
PWV SL) « L’and P(V S L) « L* occurs at about —20 dB for this 
short path. 

An obvious reason that f(—1, 0) is small for short radio links is that 
the multipath length differences are mostly less than a half-wavelength. 


VII. A DOMINANT COMPONENT INTERFERING SIGNAL 
For this case, » is bounded by 1 > wp = 1/2. By equation (36) 


awH(—1, 0) 


; L*, 1l>ped, (47) 


PV s Ll)— 
L-0 


and the corresponding power law of deep fades is 
PW sbL«L*, 1>pzh. (48) 


Since p is less than unity for this case, then as L — 0, the probability of 
deep fades decreases more slowly than square law (12). Physically, this 
means that the deep-fade problem for these links is more severe. The 
experimental data of two oversea paths (shown as curve 2 in Fig. 16.5a 
and as curve | in Fig. 16.5b of Reference 19) follow this power law (48) 
of severe fading. — 

It is true, however, that some overwater radio links still obey the 
square law (12) rather than (48).’°’” The reason is that because of the 
geometry of the radio link, the singularity in the density function 
q(R, 6) may occur, if it exists at all, at a position far away from the 
infinite fade point (R = 1, 6 = 7). Then in the neighborhood of the 
infinite fade point (R = 1, 6 = 7m), the density function q(R, 6) (or 
equivalently f(a, 8)) may still be a smooth function. 

For overland paths, it is possible that an exceptionally calm and 
stratified atmosphere would also create a stable, dominant interfering 
signal over a sustained period. Then the joint probability density 
function f(a, 8) may also be singular at the position of this stable, 
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dominant interfering signal. Therefore, the results of this section on the 
power law of severe fading may also occur on an overland radio link. 


VIII. DIVERSITY SYSTEMS 


When yz > 1, the power law of deep fades is 


PV < L) —> LO) 1, ne (49) 
PV SL)<L*, p>. (50) 


As L decreases, the probability of deep fades decreases faster than those 
following the square law (12). Physically, this means the problem of 
fading for this case is less severe than those following the square law (12). 
. The experimental data’'?"’ of composite signals of the outputs of 
diversity combining systems show that the amplitude distributions 
of composite signals in the deep-fade region obey the power-law equa- 
tion (50) rather than the square-law equation (12). 

Since » > 1 implies f(a, 8) has a zero at the deep-fade point (2 = —1, 
B = 0), these results show that the artificial active combining devices 
of diversity combining systems serve to create a zero at (a = —1, 8 = 0) 
of order 2(u — 1) of the density function f(a, 8) of the equivalent 
interfering vector of the output composite signal. The value of » depends 
on the order of diversity. By comparing the power law (50) to the 
experimental data’'’’® and the theoretical results on the diversity 
systems, we find that for most overland paths, the value of » for the 
composite signal is equal to the order of diversity. 


IX. ONE-ECHO MODEL 


In the model described by equation (15) for the fading signal, if 
there is only one echo and if the magnitude of this echo is a constant, 
then R = A is a constant rather than a random variable. For this 
idealized case, the joint probability density function q(R, 6) of R and 6 
contains a delta function 


where W(6) is the probability density function of the random relative 
phase @ between the echo Ae” and the constant vector. 

For convenience, we shall call this specialized model the one-echo 
model. (In the literature, itis also known as the two-ray model.) 

In the analysis, we shall assume that A < 1. The case where A > 1 
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can be treated similarly simply by switching the roles of echo and the 
constant vector. 

Since the joint probability density function of the interfering vector 
of the one-echo model does not belong to the class of f(a, 8) discussed 
in Section IV of this part, the results of previous sections are not 
directly applicable. Nevertheless, substituting the density function (51) 
into the general formulation (21) for P(V S L) and integrating over R 
(see Fig. 4) yields 


rtOn 
P(V<L = / W(e) ao, (52) 
w—-OL 
where 


2 2 
0, = cos” (+4 - 2) CHASES GS De 465 

Since the behavior of W(@) in the neighborhood of 6 = =z is important 
for the analysis of deep fades, we shall assume that W(@) is smooth in 
this neighborhood so that the Taylor series expansion of W(6) is appli- 
cable. Then 


WO) = We) + Wi@o- 7) + BOO t--, 64) 
where 
Wie) = FeWO| , m= 1B. 68) 


Substituting equation (54) into equation (52) and carrying out the 
integration yields 


2 2 
L*=1+2R cos 6+R* Y7R=A 





Fig. 4—The range of phase, (rt — 6,) S 6 S (x + 6,), in which V S L. 
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00 Ws(m) | = (2 as Az — ae 
PVs Ll) =2 Yoga pi Ce = pg ; 
(i+ 4)2L2(1- 4). (6) 


If the magnitudes of the two vectors are equal, then A = 1, and equa- 
tion (56) becomes 


PV SL) =2 > 4 Wasa | eos (1 _ ay 2=>L=0. 


2S + 1)! 2 
(57) 
Since 
cos’ (1 — z) = sin’ Ee — ale : (58) 
z (1 = ay Lay (59) 
IL-0 


Then the behavior of P(V < L) given by equation (57) in the deep-fade 
region is 


PY Ss L) — 2W)L; (60) 
. p(L) ad 2W(r). (61) 


This result shows that as long as the probability density function W(6) 
of the random phase @ is neither singular nor zero at @ = 7, then the 
cumulative amplitude distribution of the one-echo model with R = 
always obeys the power law 


PVsLD«L (62) 


in the deep-fade region no matter whether W(6) is uniform or not. 
Equation (57) shows that the nonuniform part of W(6) contributes 
only to the high-order terms of P(V S L) and does not affect the 
behavior of P(V S L) in the deep fade region. 

If the distribution of the random relative phase is uniform in (0, 27), 
then 


1 
WA 5, (63) 
W,(r) = 0, n=l. 
Equation (56) specialized to this case is 
ne 3 
PV sbL)= a cos" (it4—L). (64) 
T 2A 
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If A = 1, then equation (64) becomes 


P(V < L) = . cos * (1 — Z). (65) 


As far as the deep-fade region is concerned, equation (57) can also 
be written as 
PV Ss L) = 2W(i@)L + O(L'*"), n > 0; (66) 
”. p(L) = 2W(r) + O(L’), (67) 
where 0(L") is a symbol to denote the component which goes to zero 
at a rate equal to or faster than that of L” as L > 0. 

Although the one-echo model of this section and the other two cases 
discussed in Appendix C do not exhaust all the situations where H(a, 8) 
is not analytic, the main objective is to show that the assumption of 
Taylor series expansion of H(a, 8) in Section V of this part is not 
strictly necessary for the derivation of the power law of deep fades. 


To unify the representations for all the cases considered in this paper, 
we shall rewrite equations (66) and (67) as 


P(V Ss L) = 2W(a)L”* + O(L***") (68) 
and 
p(L) = 2W(@)L** + OL" **"), w=, 1>0. (69) 


Then the amplitude distributions of deep fades of all the cases discussed 
in Part 1 can be summarized as 


P(V Ss L) = d,L™ + 0(L?**") (70) 
and 
p(L) = 2nd," + 0" **"), bea > 0: (71) 


If H(a, B) is continuous at (a = —1, 6 = 0), then d, = rH(—1, 0)/u. 
If H(a, B) is discontinuous at (a = —1, 8 = 0), thend, = rH(—1, 0)/u. 
For a one-echo model, d. = 2W(m) and wp = 1/2. 
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Part 2. Expected Number of Fades and Average Fade Duration 
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I. INTRODUCTION AND SUMMARY 


In Part 1 we investigated the amplitude distribution of a fading 
signal. In a long time period, the cumulative amplitude distribution 
PV Ss L) tells us the expected fraction of this time period that the 
signal will fade below any specified signal level L. However, P(VV S L) 
does not tell us anything about the dynamic aspects of the fading 
signal. For example, a large number of short fades and a small number 
of long fades may have the same amplitude distribution. 

Some communication systems may tolerate the short fades but not 
the long fades. Furthermore, in the design of a diversity combining 
device, a distortion equalizer, or an automatic gain controlling device 
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to combat the fading problem, one needs information on the dynamic 
behavior of the fading signal. 

In Part 2 we present the results of our investigation on the expected 
number N(Z) per unit time that the signal V(é) fades below a given 
signal level L; and the average duration {(Z) of fades below L. 

The analysis is based on the general integral formulation of N(L) 
by Rice”*’” and Vigants* and our results for P(V < L) in Part 1. Again, 
we do not impose the restrictive assumption of the complex Gaussian 
model so that the theoretical results may be applied to a wide class of 
fading problems. 

In the study of N(L) and i(L), the experimental data’° for N(L) 
or i(L) are often plotted on a log scale as shown in Figs. 5 and 6. It is 
an experimental fact that the data for N(Z) and i(Z) can be well 
represented by straight lines on this kind of graph paper for fade depth 
deeper than —10 dB, as shown in Figs. 5 and 6. The slopes of these 
straight lines are directly related to the power laws of N(L) and t(L) 
in the deep-fade region. The experimental observations of N(L) and 
1(L) are summarized below: 


(z) The experimental data show that N(L) for most nondiversity 
fading signals obeys the power law, N(L) « L, in the deep-fade 
region. 

(2) For a short radio link from Villa Rica to Palmetto, Georgia, 
the N(L) of a nondiversity signal follows the cubic power law, 
N(L) « L’, in the deep-fade region. 


3 
FS 


EXPECTED NUMBER N(L) OF FADES 
2 





re) -10 -20 -30 —40 —50 
SIGNAL LEVEL £ = (20 LOG;9L) IN dB 


Fig. 5—Number of fades below signal level L. 
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102 


3 


AVERAGE FADE DURATION 
IN SECONDS 





0 -10 -20 -30 -40 -50 
SIGNAL LEVEL 2 = (20 LOGioL) IN dB 


Fig. 6—Average duration of fades below signal level L. 


(iit) The experimental data of N(L) of composite fading signals of 


most dual diversity systems also follows the cubic power law, 
N(L) « L’, in the deep-fade region. 


(iv) The available experimental data on average fade duration ¢(L) 


all obey the universal power law, t(L) « L, in the deep-fade 
region. This includes the fading signals of nondiversity systems, 
diversity systems, long radio links, and short radio links. 


In summary: 


@) 


(22) 


(ii2) 


Our theory indicates that if the joint probability density 
function f(a, 8) of the resultant interfering vector, Re’’ = a + jp, 
is a smooth function which is neither singular nor zero at the 
infinite fade point (a = —1, 8 = 0), then for small LZ 


N(L) & ra,f(—1, 0)L & L, 


where a, is a constant approximately equal to the average 
positive derivative of the amplitude of the fading signal in the 
deep-fade region. 

If f(a, 8) is singular or zero at (2 = —1, B = 0) then 


N(L) © wH(-1, 0)a,L7*"? «© L?*7 


for » = 1/2 and small L. The cubic power law, N(L) « L’, of 
dual diversity systems and short radio links can be explained 
by this result when » = 2, which means f(—1, 0) is zero or is 
negligibly small. 

The theory predicts that the average fade duration always obeys 
the power law ¢(L) « L for small ZL no matter whether f(a, 8) is 
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smooth, singular or zero at (a = —1, 8 = 0). This means the 
power law, i(Z) « L, is invariant with respect to variations of 
fading environment and diversity combinations of fading signals. 
This prediction agrees with the available experimental data. 

(tv) The theoretical results on N(Z) and i(L) for the one-echo model 
are shown to be incompatible with the experimental data of 
most overland microwave radio links. Therefore, the one-echo 
model is not suitable for the study of the statistics of fading 
signals of these radio links. 

(v) For line-of-sight radio links at 4 GHz and 6 GHz, the average 
positive derivative of the amplitude of the fading signal is esti- 
mated to range from 2 X 107° to4 X 107° times V,,.; per second, 
where V..¢ 18 the signal level when there is no interference. 


II. GENERAL FORMULATION FOR NUMBER OF FADES 


The general expression for the expected number of fades per unit time 
of a random signal V(é), below signal level L has been shown”’’”’* to be 


NC) = f.” VoCV, ¥) Ira dV, 72) 


where V = dV/dt, and p(V, V) is the joint probability density of 
V and V. For the sake of completeness, a brief derivation of (72) is 
included in Appendix E. 

The joint probability density function p(V, V) can be written in 
terms of conditional probability” as 


pV, V) = pi(V | V)p.(V), (73) 


where p,(V | V) is the conditional probability density of V under the 
condition that the signal level is V; and p.(V) is the probability density 
of V. Substituting (73) into (72) yields 


NQ) = pt) [Voc | Dav. (74) 


Let us define 


V.(L) = 2 . ~ Vo.(V | L) av. (75) 


The physical meaning of the definition (75) is that V,(L) is the condi- 
tional average positive derivative of V under the condition V = L. 
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The factor 2 in (75) is based upon the assumption that p,(V | L) is 
symmetric about V = 0. 

Since V(¢) is a random fading signal, at a given signal level L, the 
value of its time derivative V is also random. In general, the conditional 
average positive derivative V,(L) is a function of signal level L. 

By using definition (75), equation (74) becomes 


N(L) = 3V.(L)p2(L). (76) 


Or, equivalently, 
ND) = 3.) PW SD, (77) 


which indicates N(L) proportional to the conditional average positive 
derivative of the fading signal and to the probability density of fades 
at V = L. 

II. GENERAL FORMULATION FOR AVERAGE FADE DURATION 


In a long time interval, T, containing a large number of fades,* the 
expected total length of time that the random signal V(t) spends below 
a specified signal level L is 


((L) = TP(V < L). (78) 


The expected number of fades below L.in this interval T is TN(L). 
Therefore, the average duration of fades below L is 


(L) _ PV SL) 





i) = NG) ~~ ND) (79) 
Substituting (77 ) into (79) yiel ds 
ur) = 1 PP SD i 


VAD) Pv gD) 


Equation (80) shows that the average fade duration is inversely pro- 
portional to the conditional average positive derivative of the fading 
signal at V = L. 


*In our experiment on line-of-sight radio links, the typical time interval 7 is 
a whole summer of more than 100 days in which there are more than 500 fades 
below —10 dB relative to the nonfaded signal level. 


FADING SIGNAL STATISTICS 3235 


IV. ASSUMPTION ON CONDITIONAL AVERAGE POSITIVE DERIVATIVE 


Most existing theoretical work on N(L) assumes that V and V are 
independent so that V,(L) becomes a constant which is independent of 
signal level L. In this paper we include the situation for which V and V 
are dependent and assume that V,(L) can be expanded into a Taylor 
series in the deep-fade region (i.e., small ZL); then 


where 


V.~L)=at+aLh+al? + aL? +---, (81) 
ay = lim ¥.(L) = ¥.(°), 2) 
a, = jim 3r = V.(D), (83) 
a, = 5 Hi lim 4% V.(L), etc. (84) 


The justification for this assumption is not trivial, and includes the 
following considerations: 


(7) 
(2) 


The theoretical results based on this assumption agree with the 
available experimental data. 

For a complex Gaussian model, the conditional probability 
density function p,(V | L) is known. Then with the help of the 
work of Rice”*’” the integration indicated in (75) for V,(L) can 
be carried out in closed form. These explicit expressions are 
discussed in Section III of Part 3. The results of this model show 
that if the power spectrum of the Gaussian noise is symmetric 
with respect to the frequency of the sine wave (i.e., the constant 
unit vector), then V,(L) = a, is a constant independent of L. 
On the other hand, if the power spectrum of the Gaussian noise 
is asymmetric with respect to the signal frequency, then V,(Z) 
is a function of Z and the nonconstant terms in equation (81) 
cannot be omitted. 

The theoretical work of Clarke,” Ossanna,*° and Gans® on 
mobile radio indicate that the power spectrum of the fading 
signal is generally asymmetric with respect to the received 
carrier frequency unless the straight line joining the base station 
and the mobile antenna is perpendicular to the velocity of the 
mobile and the antenna pattern is symmetric with respect to 
this line. Therefore, the work on asymmetric power spectrum, 
and hence nonconstant V.(L), is not purely academic. 
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(wii) It is known” that the correlation between any real random 
variable e(t) at instant ¢ and its time derivative ée(t + £) at 
instant (t + £) vanishes if ¢ = 0. This is often used to support 
the assumption that V and V are independent and hence V’,(L) 
is a constant. However, we know that the vanishing of correla- 
tion between V(t) and V(t + ¢£) at ¢ = 0 does not imply the 
independency of V and V unless V is normally distributed. 
Including high-order terms in equation (81) removes the assump- 
tion of independency of V and V and enlarges the applicable 
scope of this theory. 

(iv) In equations (82), (83), and (84) we define those coefficients of 
the Taylor series as the limits of ¥,(L) and its L-derivatives 
at L = 0* from the positive side. The reason is that V(é) is the 
absolute value of a fluctuating complex signal; i.e., 


V® =|e+iy| = |Ve*]. 


When the complex fluctuating signal V(é)e’*"” crosses zero, its 


absolute value V(t) may have a cusp at V = 0 as shown in Fig. 7. 
Therefore, the derivatives of V(t) may not be well defined 
at V = 0. However, the limits of the derivatives at V = 0* 
from the positive side are well defined. 

(v) Fig. 7 also shows that although V = 0 is a minimum of V(), 
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Fig. 7—The amplitude V(t) may have a cusp at V = 0 even though the complex 
fading signal is a smooth time function. For convenience, the complex fading signal 
is plotted as a real function in this example. 


FADING SIGNAL STATISTICS 3237 


V(t) = d/dt V(é) may not be zero at V = 0°. Therefore, in 
equation (81), the constant term a, = V,(0*) does not vanish 
for most cases, and cannot be omitted from equation (81). 


V. POWER SERIES FORM OF N(L) 


Substituting the power series (35) for (0/9L) P(V S L) of Part 1 and 
the power series (81) for V,(Z) into equation (77) yields 


N(L) = rH(—1, 0)aL** + eH(—1, 0)aL7* + --- (85) 


In the deep-fade region, the leading term dominates the power series 
(85). Therefore, 


N(L) & rH(-1, O)a.L7"", forsmall LZ; (86) 

. N(L) « L**, forsmall L. (87) 

Equations (85) and (86) show that if a, # 0 and if » < 1/2, then 
N(L) > « ae L — 0. This seems to be nonphysical. Therefore, we 


require that 1 = 1/2. This is the reason we impose this condition on the 
order of singularity of f(a, 8) in Section IV of Part 1.* 


5.1 Prevailing Power Law of Number of Deep Fades 
For the nondiversity fading signals of most radio links, the probability 
density function f(a, 8) of the resultant interfering signal is a smooth 
function which is neither singular nor zero at (2 = ~—1, 8 = 0). Then 
= 1 and f(a, 8) = H(a, 8). Equations (85), (86), and (87) under this 


condition become 


N(L) = rf(—1, O)a.L + rf(—1, O)a,L? 


+ [2a) d,s + rf(—1, O)a,]L? + --- (85’) 
N(L) = rf(—1, 0)aoZ, forsmall L. (86’) 
N(L) « L, for small L. (87’) 


It is seen that as long as » > f(—1, 0) > 0, then the expected number 
N(L) of deep fades always obeys the prevailing power law (86’). With 
reference to Fig. 5, the straight lines corresponding to the power law (87’) 
have inverse slopes of 20 dB per decade. 


* Tf we assume that ay = 0, then the only constraint on » is p > 0 due to the unity 
total probability. pee assuming @ = 0 implies that the time derivative of V 
is always zero at V = 0. Such an assumption is unreasonable for multipath inter- 
ference fading, but, may be useful if the effects of random circuit interruptions, such 
as equipment failure, are included in the signal fading problem. 
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Apparently, the simple condition is easily satisfied by the nondiversity 
fading signals of most microwave radio links because the experimental’~° 
data of N(Z) are mostly characterized by the inverse slope of 20 dB 
per decade in the deep-fade region. 

As L increases from zero towards unity, equation (85’) indicates 
that there may be a transition point beyond which high-order terms 
become significant and the slope begins to deviate. 

For most microwave radio links, the transition points of N(L) seem 
to be well above —10 dB. However, our latest experimental data show 
that the first transition point of N(L) of a relatively short path (path 
length 15.87 miles, f = 4 GHz) is below —20 dB. The inverse slope 
of N(L) in the region from — 10 dB to —20 dB is approximately 20/3 dB 
per decade of number of fades. This indicates that the third term 
[2a, ds + masf(—1, 0)JL* dominates in the region 0.8 > L > 0.1 for 
this path. 


VI. AVERAGE DURATION OF DEEP FADES 
Substituting the power series (31) of P(V S L) of Part 1 and the 
power series (85) of N(Z) into equation (79) yields 


aH(—1, 0) 5 Ps gd Te"? aa + . 


earner 9 aera -7 aa FO (88) 
rH (—1, 0)a,L7""! + rH(—1, 0)a,L7* + 


In the deep-fade region, equation (88) becomes 


iL) = = L, forsmall L; (89) 
Mao 
. UL) «< L, forsmall JL. (90) 


On Fig. 6, the straight lines corresponding to the power law (90) have 
an inverse slope of 20 dB per decade of fade duration. The experimental 
data’~° agree with this conclusion on the slope of é(L) when plotted 
on Fig. 6. 


VII. INVARIANCE OF POWER LAW OF AVERAGE FADE DURATION* 


Equations (36), (87), and (90) show that in general P(V S L), N(L), 
and i(L) obey the following set of power laws of deep fades: 


_ *In equation (26) of Ref. 32, Rice has already predicted that the power law 
i(L) « L for small L may be applicable to cases more general than the complex 
Gaussian model even though most of his work in Ref. 32 is devoted to the statistics 
of a sine wave plus a narrowband Gaussian noise. 
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P(V SL) « L* | (91) 
N(L) « Le for » 24 andsmall L. (92) 
iL) « L (93) 


It is seen that the power laws of P(V < L) and N(L) depend on the 
value of » which depends on whether f(a, 8) is smooth, singular, or zero 
at (a = —1, B = 0). On the other hand, the power law (93) for the 
average fade duration is invariant with respect to yp. Since the behavior 
of f(a, 8), and hence the value of u, depends on fading environment, 
we conclude that the power law (93) for the average fade duration is 
insensitive to the fading environment in contrast to the power laws 
of P(VV < L) and N(L). 

In equation (89), notice that ua, , and hence ¢(L), does depend on the 
fading environment. However, it is the power law, i(L) « L, which is 
insensitive to the fading environment. 

Vigants,*’°’"* Crawford, Hogg, and Kummer” have investigated the 
effects of diversity on P(V < L), N(L) and i(L). The theoretical results 
and the experimental results of these authors show that in the deep-fade 
region, the diversity drastically changes the power laws of P(V S L) 
and N(L) but does not affect the power law i(L) « L. For example, the 
results of Vigants are shown in Table I. From this table it is seen that 
the power laws of P(V S L) and N(ZL) depend on the diversity combina- 
tion of fading signals, whereas the power law, i(L) « L, of average fade 
duration is invariant. 


VIII. INCOMPATIBILITY BETWEEN ONE-ECHO MODEL AND OVERLAND RADIO 
LINKS 


_ Equations (60), (61), (91), (92), and (93) show that p(L), N(Z), and 
t(L) of the one-echo model (with equal magnitudes, A = 1) in the 


TaBLE J—Errects oF DIversiry oN PowreR LAws 
oF DEEP FADES 


Nondiversity Diversity 
PV SL) i? (1/q) L4 
N(L) cL (2c/q) L’ 
i(L) (1/e) L ies L 


Remark: In this table, the parameter, c, as defined by Vigants, is equal to ao 
of this paper; and the parameter, g, as defined by Vigants, is equal to 2/rH(—1, 0) 
of this paper. 
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deep-fade region are 


p(L) « al (94) 
N(L) « L’¢ forsmall L. (95) 
i(L) « a (96) 


On the other hand, the long-term experimental data of a nondiversity 
signal of most overland microwave radio links indicate that 


p(L) « L (97) 
N(L) « Le forsmall L. (98) 
i(L) « ‘ (99) 


The experimental results (97) and (98) disagree with (94) and (95) of 
the one-echo model. 

In view of this disagreement, we may want to check the effect of the 
assumption (81) of V.(Z) on the theoretical results of the one-echo 
model. Although we know that the constant term a, of (78) generally 
does not vanish, yet we may deliberately set a, = 0 and see what kind 
of theoretical results we get. 

If we do so, the theoretical results of the one-echo model become 


p(L) « L (100) 
N(L) « L forsmall JL. (101) 
i(L) « =| (102) 


Under this modified assumption, p(L) and ¢(L) of the one-echo model 
disagree with the experimental results (97) and (99). Similarly, forcing 
the coefficients of other higher order terms of (81) to zero also yields 
theoretical results which disagree with the experimental results. There- 
fore, we conclude that the one-echo model is not suitable for the study 
of the fading signals of most overland microwave radio links. 

However, we emphasize that the experimental data mentioned in this 
section are restricted to the long-term data of overland microwave radio 
links. Therefore, the incompatibility of the one-echo model with these 
data does not necessarily exclude the use of this model for the study of 
other fading problems. 
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TX. APPROXIMATE AVERAGE POSITIVE DERIVATIVE V5 OF LINE-OF-SIGHT 
MICROWAVE RADIO LINKS 


For line-of-sight microwave radio links, let V,., be the signal level 
when there is no interference, and V;,a(t) be the random fading signal 
when the interference appears. In our analysis 


VQ = Leal (103) 
Veet 
is a normalized fading signal. : 

By comparing the experimental data of ¢(L) and the theoretical 
equation (89) for ¢(L), we can estimate the value of a, of the radio link. 
The value of » in equation (89) can be determined from the experimental 
data on the power laws (36) and (86) for P(V < L) and N(ZL) of the 
same radio link. 

Our experimental data of several line-of-sight microwave radio links 
in Ohio and Georgia indicate that the value of a, ranges from 2 X 10° 
to 4 X 10°. In the deep-fade region where ZL is small, equation (81) 
shows that V,(Z) = a,. 


~a=2X10°~4X 10°. — (104) 


average 








d 
qi VO 
Substituting (103) into (104) yields 


~(2X 10° ~4 X 10°): Veer . (105) 


average 


d 
dt V raat) 








Thus, the average positive derivative of the unnormalized fading signal, 
V.a(t), of these microwave radio links ranges from 2 X 107° to 4 X 10° 
times V,., per second. 

These approximate values of average positive derivative are valid 
only in the deep-fade region because they are deduced from the experi- 
mental data of deep fades. The path length of these radio links ranges 
from 15 miles to 36 miles. The operating frequencies are in 4-GHz 
and 6-GHz bands. 


X. GENERALIZED ASSUMPTION ON Y,(L) 


In Part 1, we indicated that the assumption of the Taylor series 
expansion of H(a, 8) is not strictly necessary for the validity of the 
power law of P(V S L) deduced from the experimental data. In this 
section, we point out that the assumption of the Taylor series expansion 
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of V,(L) in Section IV of this part is also not strictly necessary for the 
validity of the power laws of N(L) and ¢(L). From a theoretical view- 
point, the assumption of Y,(Z) can be generalized to the following form: 


V.(L) =a, +0(L"), >0 (106) 


where 0(L") is a symbol to denote the component which goes to zero 
at a rate equal to or faster than that of L’ as L > 0. 


In assumption (106), we do not require the existence of the limits 
in equations (83), (84), etc. Therefore, the assumption (106) is less 
restrictive than the assumption (81). It can be shown that the power 
laws of N(L) and ¢(L) of deep fades based on (106) are the same as 
those based on (81). However, at the present time, we do not have any 
practical evidence to necessitate the use of (106). Therefore, we merely 
point out the possibility but do not explicitly carry out this generalized 
analysis. 
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Part 3. Special Topics on Statistics of a Fading Signal 
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I. INTRODUCTION AND SUMMARY 


In Parts 1 and 2, the analysis was oriented towards an explanation 
of the experimentally observed common behavior of a fading signal Ve’*. 
The basic assumptions of the theoretical model are kept to a minimum 
in order to include the widest possible variation in practical fading 
environments. During the development of this general analysis, we have 
gained a new insight into several topics related to fading signals as 
investigated by previous authors. 
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Part 3 of this paper is a collection of theoretical treatments of several 
special topics relating our generalized analysis to the work of previous 
authors. These topics include the sum of n unit vectors with random 
phases; a sine wave plus Gaussian noise, m-distributions, chi-distribu- 
tion, Rayleigh distribution, and log normal distribution. 


In summary: 


(2) 


(72) 


(ii2) 


(wv) 


In Section IJ of this part, the results of Part 1 are applied to 
find the amplitude distribution of the sum of n unit vectors with 
uniformly distributed random phases. For n = 3, the analysis 
shows that the amplitude distribution always follows the square 
law P(V < L) « L’ for small L. On the other hand, when 
n = 2, the amplitude distribution follows the power law 
PV s L) « L for small LD. 

In Section III, we investigate the model of a sine wave plus a 
narrowband Gaussian noise for the fading signal. By using the 
closed-form solutions of Rice, it is shown that if the power 
spectrum of the Gaussian noise is symmetric with respect to 
the frequency of the sine wave, then the amplitude V and its 
time derivative V are independent; and the conditional average 
positive derivative V,(L) is a constant. On the other hand, if 
the power spectrum is not symmetric, then V and V are de- 
pendent; and the conditional average positive derivative V,(L) 
is a function of signal level V = L. 

As an example, the fading signal spectral density of a mobile 
radio is generally not symmetric with respect to the received 
carrier frequency. Therefore, in the analysis of N(L) and ¢(L), 
it is not safe to assume that V and V are always independent. 
In Section IV, we investigate the theoretical condition (147) 
on the joint probability density function f(a, 8) of the interfering 
vector such that the amplitude distribution of the fading signal 
belongs to the family of m-distributions which includes normal 
distribution, Rayleigh distribution, Maxwell distribution, and 
all of chi-distributions as special cases. 

It is also shown that the set of m-distributions behave like a 
log normal distribution within a small range (148) of signal 
level near its rms value. This result shows that in the interpre- 
tation of the experimental data, one must be cautious in attempt- 
ing to estimate the tails of the distribution by an extension from 
the middle section of the distribution. 

We find that 1 in general, the integral transformation (20) from 
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f(a, 8) into P(VV S L) is not unique. Physically this means the 
signals of fading environments with different f(a, 8) can have 
the same amplitude distribution P(V S L). As an example, this 
nonuniqueness shows that specifying a Rayleigh distribution 
for P(V Ss L) does not necessarily imply that there are a large 
number of interfering signals; nor does it necessarily imply that 
the real and imaginary parts of the fading signal are normally 
distributed with zero mean. 


II. SUM OF 7” UNIT VECTORS WITH UNIFORM RANDOM PHASES 


The amplitude distribution of the sum of two unit vectors with 
uniformly distributed random relative phase has been shown in Sec- 
tion IX of Part 1 to be, forO = L S 2, 


PWV <L) 


2 
: cos * (1 — Z) : (107) 
T 2 


and 
1 


In the deep-fade region where L is small, this amplitude distribution 
obeys the power law 


(108) 


ale 


pL) = PV SD) = 


PV SL) « L. (109) 


The sum of n unit vectors with uniformly distributed random phases 
has been investigated previously by many authors.**'**~*” The mathe- 
matics involved in obtaining the amplitude distribution for any arbitrary 
n = 3 is fairly complicated. Computer numerical integration is needed 
to show the distribution explicitly. In this section, we shall avoid the 
complicated mathematics and shall apply the results of Part 1 to show 
that the amplitude distribution for any arbitrary n = 3 in the deep-fade 
region always follows the square law: 


P(V SL) « L’, for small L. (110) 


The sum of n unit vectors with random phases can be written as 


Ve" = eC (111) 
i=1 


E +. Fen | 79, (112) 
1=2 
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= [1 + Re*Je’", (113) 
where 
Re’ = Se gf tint, (114) 
i=2 


Rosenbaum* has indicated that if all the phases {6,;}i=? of the unit 
vectors are independently and uniformly distributed in (0, 27), then 
Ve'* and Re’’ have circular symmetric probability density functions; i.e., 
the amplitude and the phase are independent and the phase is uniformly 
distributed in (0, 27). 

It then follows that the random signal represented by equation (113) 
for any arbitrary n 2 3 is a special case of Appendix B. The case for 
nm = 2 is an exception because the joint probability density function 
q(R, 6) contains a delta function whereas the q(R, 6) in Appendix B is 
assumed to be a smooth function. 

The sum of three unit vectors can be considered as a unit vector 
suffering interference by a random vector R(t)e”’“” which is the sum of 
the other two unit vectors. The amplitude distribution g(R) of R(t) is 
given by equation (108) except for the replacement of the notation L 
by R. Equation (108) implies 


g(l) = (115) 


Therefore, g(R) for this case is a smooth function which is neither 
singular nor zero at R = 1. Then equation (165) shows that 


ml a a. 2 2 
PV SL) =i3g0)L = aa FSi aa Bae (116) 
The sum of n unit vectors, n = 3, can be considered as a unit vector 
suffering interference for a random vector R(t)e’’“” which is the sum of 
the other (n — 1) unit vectors. It is obvious that g(1) # 0 simply 
because each of the (n — 1) unit vectors has unity amplitude. Then, 
the results (see Appendix B) imply 


PWV < L) « L’, forsmall Land n = 3. (117) 


On a log-versus-dB graph paper, as shown in Fig. 1, the power law (109) 
implies a straight line with an inverse slope of 20 dB per decade of 
probability whereas the square law (117) implies a straight line with 
the same inverse slope of 10 dB per decade of probability as that of 
Rayleigh distribution. 

Therefore, we conclude that for n = 2, the distribution of deep fades 
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is characterized by the inverse slope of 20 dB per decade of probability 
whereas for any n = 3, the distribution of deep fades is always charac- 
terized by the inverse slope of 10 dB per decade of probability. This 
conclusion agrees with the numerical results of Norton, et al., in Fig. 2 
of Reference 33. 


Ill. A SINE WAVE PLUS A GAUSSIAN NOISE 


The statistical behavior of a sine wave plus narrowband Gaussian 
random noise has been investigated in great detail by Rice.” ”**? In this 
section we shall apply our general analysis to this case to show the 
consistency of our results with the work of Rice. Furthermore, we shall 
also use the closed-form solution of N(Z) and t(Z) obtained by Rice 
to show that the conditional average positive derivative V..(Z) can be 
either a constant or a function of signal level LZ, depending on whether 
the power spectrum of the noise is symmetric or asymmetric with 
respect to the frequency of the sine wave. 

In this model, the interfering vector, Re’? = a + j8, represents the 
envelope of a narrowband Gaussian noise; the constant vector represents 
the sine wave with a constant amplitude and frequency f,. The joint 
probability density function f(a, 8) is a two-dimensional normal density 
function; Le., 


fla, B) = 525 exp [—(@ + 6°)/2b), (118) 


where a and 6 are assumed to be independent normal random variables 
with the same variance b, and zero mean. 
The well known Rice distribution for the amplitude of this model is 


pL) = 2 PW < l= 7 (22) ‘sie (==) , (119) 





where I,(~) is the modified Bessel function of zero‘ order, and Q is 
the magnitude of the sine wave. In our analysis, Q = 1 because all the 
signals are normalized to the magnitude of the constant vector. 

Rice” has also shown that the joint probability density function 
p(V, V) for this model is 


pV, V) = p {sat [B(V? — 2VQ cos + Q’) 


Pas [ex 


+ (boV + 6,Qsin a} dp, (120) 
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where ¢ is the phase of the resultant fading signal Ve’*; f, is the frequency 
of the sine wave, w(f) is the power spectrum of the Gaussian noise, and 


b= ny [wt af, 2 =0,1,2, (2 


b= bobs = bi . (122) 
3.1 Amplitude Distribution of Deep Fades 


The normal density function (118) is obviously a smooth function 
which is neither singular nor zero at the infinite fade point (a = —1, 
8 = 0). Then the results of Section VI of Part 1 predict that the ampli- 
tude distribution of the fading signal in the deep-fade region is 


PWV < LD) &crf(-1,0L’ (123) 
= 5, exp: (3) x. (124) 


On the other hand, the limiting form of the Rice distribution for small L 
is 


p(t) = © exp (SX) ; (125) 


a pendeets 3f).r. for L«KQ=1. (126) 
2b 2b 

It is seen that our result (124) agrees with the Rice distribution in the 
deep-fade region. 

The square law (124) implies that on a log-versus-dB graph paper, 
the Rice distribution in the deep-fade region is always characterized 
by the prevailing inverse slope of 10 dB per decade of probability. The 
numerical results of Norton, et al., in Fig. 5 of Reference 33 agree 
with this prediction. 


3.2 Symmetric Power Spectrum and Constant V.(L) 


If the power spectrum w(f) of the Gaussian noise is symmetric about 
f,, then b, = O and the integration of (120) under this condition yields 


10. = {os [EE wo [==] 0 


= pV | V)p.(V) (128) 
= p,(V)p.(V), (129) 
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where 
es (24) |= yr =) 
p(V) = i I, by ) oxP 3b, (130) 
is the Rice distribution for V; 
| Sr] 
V 131 
pV) = Fae exw | (131) 
is a normal density function for V; and 
B 5 
bb = be = Vins When 0b, = 0. (132) 
0 


Equations (127) and (129) show that if the power spectrum is symmetric 
about f, , then the envelope V and its time derivative V are independent, 
and V is normally distributed. 7 
Substituting (131) and (132) into the definition (75) for V,(Z) yields 
rn QV 
V.(L) = —' = a = constant. 133 
+(Z) Soe 0 (133) 
Therefore, the conditional average positive derivative 7,(L) for this 
model is a constant if the power spectrum is symmetric about f, . 
Substituting (180) and (133) into the general expression (76) for 
N(L) yields 








N(L) = Se 7 (22) Sm |= aaa Q | (134) 
= $aop2(L). (135) 
Then 
2 @(V) dV 
PV SL) _ af . 
A= Dy pay ue 
In the deep-fade region where L K Q = 1, 
N(L) & woe Gx * exp | 5 | if (137) 
= ra f(—1, 0)L, (138) 
WL) & +L. (139) 


Qo 
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It is seen that equations (138) and (139) agree with equations (86’) 
and (89) of Part 2. 


3.3 Asymmetric Power Spectrum and Nonconstant V.(L) 


If the power spectrum w(f) is not symmetric about f, , then V and V 
are dependent and b; # 0. The joint probability density function 
p(V, V) for this case cannot be written as the product of the individual 
probability density functions of V and V. For this case, Rice” has 
obtained N(L) by substituting (120) into the general expression (72) 
and carrying out the integration. This gives 





B 
p(L) do n=0 1 (= ber’ 





NE SS Lo 
_ 2 oni \ QL 
anf) BP 
on) bor” (2) 
| 1(2 + QL © | ae ’ (140) 
where I,,(-~~) is the modified Bessel function of order n, and 
2 _ b,Q’ 
OS pe (141) 


Comparing equation (140) and the general oapieeon (76) for N(L) 
shows that 


V.(L) = 








pee ee ao 
= eee 
bo 


[HQ (2) 


It is seen that when V and V are dependent, then the conditional 
average positive derivative V.(L) is a function of signal level V = L. 

The expected number of fades N(L) and the average fade duration 
t(L) for this case in the deep-fade region are 


N(L) & raof(—1, 0)L 
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where 


a = lim VAL) = (28) 3 ee E x 2(n ry me 


The work of Clarke,” Ossanna,”” ‘and Gans*' on mobile radio indicates 
that the power spectrum of the fading signal is generally not symmetric 
with respect to the received carrier frequency unless the straight line 
joining the base station and the mobile antenna is perpendicular to the 
velocity of the mobile and the antenna pattern is symmetric with 
respect to this straight line. Therefore, in the theoretical work of N(L) 
and i(L), it is not safe to assume that V and V are always independent. 


IV. M-DISTRIBUTIONS, CHI-DISTRIBUTIONS, AND RAYLEIGH DISTRIBUTION 


In the study of the experimental data of amplitude distributions of 
short-term high-frequency long-distance propagations, Nakagami*® 
found that the set of experimental data can well be described by a 
family of m-distributions:* 


p(L) = om aa — 12" exp |= mle | (143) 








where Q is the mean square value of the fading signal. The operating 
frequency ranged from 10 MHz to 20 MHz and the path length from 
1500 kilometers to 9000 kilometers. Nakagami indicated that these 
results were obtained from short records of data from three to seven 
minutes in length in order to avoid the effects of slow fading on the 
distribution of rapid fading. 

The various properties of the m- -distributions have been investigated 
in detail by Nakagami.*” It is easily shown that the set of chi-distribu- 
tions” is a subset of m-distributions by setting 2m = any positive 
integer in (143). This means the normal distribution, Rayleigh dis- 
tribution, and Maxwell distribution are also special cases of m-distri- 
butions when m = 1/2, 1, and 3/2 respectively. On Rayleigh paper, 
all the m-distributions appear to be straight lines passing through the 
common point of 50 percent at 0 dB, with different slopes which depend 
on the value of m. The graphical representation of m-distributions can 
be seen in Reference 39. 

However, in Reference 39, one does not know the theoretical condi- 
tion under which the amplitude distribution of a fading signal will 


*'To avoid possible confusion, we emphasize that the Nakagami distribution 
mentioned in Section I of Part 1 is not the m-distribution discussed in this section. 
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belong to this family of m-distributions. In this section we shall find 
the condition on the joint probability density function f(a, 6) of the 
interfering vector such that the amplitude distribution P(V < L) will 
belong to m-distributions. 

Expanding the exponential function in equation (143) into a power 
series gives 


o) = Tahar 2 Pa ary aa oe 


Comparing equation (144) and the general power series (35) for p(L) 
shows that 


B= Mm”, (145) 
: 2(—1)* 
OS + 2) dress = TO (146) 
(2) T(u)-S! 
Bb 
Substituting (84) into (146) gives 

v=8 = Ss 

Tv Hos-2y,2v(—1, 0) = ( 1) S= 0, 1, 2. state (147) 


28 £4 (VIS — »)! “P20 


Thus, equations (145) and (147) are the general conditions on the 
interfering vector such that P(V < L) is an m-distribution. 

In Section V of Part 2 we showed that » 2 1/2 whichimplies m 2 1/2. 
Nakagami*® has also found this condition on the parameter m by a 
different approach. Since m can be any value 2 1/2, equation (143) 
represents an infinitely large family of distributions. 


4.1 Log Normal Behavior of m-Distribution Near the RMS Value 


The experimental data of optical propagation*® *” and line-of-sight 


radio links show that the distributions of the signal scintillation near 
its average value are approximately log normal. Usually the accuracy 
of experimental data is best in the middle section of the distribution 
and deteriorates towards the tails. It is quite tempting to estimate the 
tails of the distribution by an extension from the middle section. The 
deviations of the experimental data at the tails are often attributed to 
the experimental error. 

However, de Wolf*® and Deltz and Wright’” have pointed out that 
the use of the middle section of a log normal paper may not be a reliable 
test of the log normal distribution. The differentiation between the log 
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normal distribution and certain other distributions may be significant 
only at the tails rather than the middle section of these distributions. 

Nakagami’’ has pointed out that all the m-distributions behave like 
a log normal distribution for the fading signal V (¢) in the neighborhood 
of its rms value cg = VQ. The explicit bounds on the signal level 
within which this approximation holds are: 


20 





L 
It is shown in Appendix F that the m-distribution within the signal 
range (148) is approximately equal to 


p(L) = Een | exp [—2m(In L — In V0)"), (149) 


which is a log normal distribution for the signal level L. 

Therefore, the m-distributions, including the normal distribution, 
Rayleigh distribution, Maxwell distribution, and chi-distributions, all 
behave like a log normal distribution within the signal range (148). This 
result points out that in the interpretation of the experimental data, 
one must examine the behavior of the data not only inside but also 
outside of the range (148) in order to assert their distribution. 


V. NONUNIQUE RELATION BETWEEN AMPLITUDE DISTRIBUTION AND f(a, 8) 


From the general integral relation (20) between P(V S JL) and 
f(a, B), it is seen that any component of f(a, 6) that is antisymmetric 
with respect to (1 + a) and/or 6 will cancel out in the integration (20), 
and contributes nothing to P(V < L). This means that there are many 
different f(a, B)’s, with the same symmetric* part and different anti- 
symmetric* part, which correspond to the same amplitude distribution 
PV s L). 

Mathematically, this means the integral transformation (20) from 
f(a, B) into PV S L) is not unique. Physically, this means the fading 
signals in fading environments with different f(a, 8) can have the same 
amplitude distribution. 

Furthermore, even if we restrict f(a, 8) to functions symmetric with 
respect to (1 + a) and 8, the relation between f(a, 8) and P(V S L) is 
still not unique. We shall demonstrate this nonunique relation specifi- 
cally by using the results for the m-distributions previously discussed. 
We notice that for each S, equation (147) is an algebraic equation 


* With respect to (1 + @) and £. 
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for (S + 1) unknowns {Hos-2,,2,(—1, 0)}"2§ . It is then obvious that. 
there are infinitely many different sets of {H2s-2,.2,(—1, 0)}"2% which 
will satisfy equation (147) because there is only one equation for 
GS + 1) unknowns. This nonuniqueness gives a great freedom for the 
wide variations of the individual term, Hes-2,,2,(—1, 0), which is the 
even-order partial derivative of H(a, B) at (a = —1, B = 0). From the 
Taylor series (27), it is seen that Hes-2, 2,(—1, 0) is the coefficient of 
the even-order term (1 + a)?*~*’8”” which is symmetric with respect 
to (1 + a) and g. Therefore, the relation between f(a, 8) and P(V S L) 
is not unique even if f(a, 8) is symmetric with respect to (1 + a) and 8. 

A more detailed discussion of this nonunique relation in polar co- 
ordinates is given in Appendix G. 


VI. PHYSICAL MODEL AND RAYLEIGH DISTRIBUTION 


In this section, we shall show that specifying a Rayleigh distribution 
for the amplitude of a complex fading signal Ve’* does not necessarily 
imply that there is a large number of interfering signals; nor does it 
necessarily imply that the real part and the imaginary part of the 
fading signal are normal with zero mean. 

Let x and y be the real part and the imaginary part respectively of 
the complex fading signal Ve’*, and let F(z, y) be the joint probability 
density function of x and y. Since V? = x + y’, then the probability 
of V s Lis the probability of x and y falling within the circular region 


gt+ty <L’. (150) 


Therefore, P(V S L) is the integration of F(z, y) over the circular 
region (150); i.e., 





y=L zaVEI—y? 
P(V <L) = i i F(z, y) dx dy. (151) 


aa retry 
6.1 Number of Interfering Signals 


In the most common derivation of the Rayleigh distribution, the 
fading signal is assumed to consist of a large number of random, inde- 
pendent interfering signals, 


Ve® = DoE =xt+ yy. (152) 
i=1 


Furthermore, it is assumed that none of the components {#;}!=7 pre- 
dominates in the summation (152). Then by the central limit theorem, 
one argues that as the number, , of interfering signals approaches 
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infinity, the real part x and the imaginary part y of the fading signal 
Ve'* become independent normal random variables with the same 
variance and zero mean. This implies that V and ¢ are independent and 
¢ is uniformly distributed in (0, 27). Under this condition, the distribu- 
tion of the random amplitude V is Rayleigh. 

However, by an observation similar to those in Section V of Part 3 
and Appendix G, we realize that the transformation (151) from F(a, y) 
into P(V S L) is not unique. Given an amplitude distribution P(V S L), 
there correspond infinitely many different F(z, y)’s. In other words, 
the independent normal distribution for x and y with the same variance 
and zero mean is only a sufficient condition but is not a necessary 
condition for the amplitude V to have a Rayleigh distribution. Since 
xz and y do not have to be normal, then the number n of the interfering 
signals does not have to be large. 

In Appendix G, we have shown that the relation between P(V S L) 
and F(x, y) becomes unique if the following two additional conditions 
are imposed; 


(<) V and ¢ are independent, and 
(it) $ is uniformly distributed in (0, 27). 


For long radio links such as beyond-the-horizon radio links, the condi- 
tions 7 and 2 seem to be applicable. However, for line-of-sight radio 
links, our experience indicate that the phase ¢ has much higher tendency 
of wide variation during the deep fade where V is small. This means 
for short radio links, V and ¢@ may not be independent and ¢ may not 
be uniform. Therefore, in our general analysis we do not impose the 
conditions 7 and 77. 


6.2 Mean Values of x and y 


In the integral relation (151) the antisymmetric part of F(z, y) 
contributes nothing to the amplitude distribution P(V S LD), but does 
affect the mean value of x and y. Then by adding a suitable* antisym- 
metric function to F(x, y), the mean values of x and y can be changed 
arbitrarily without affecting the amplitude distribution P(V < L) of 
the fading signal Ve7®. 

In other words, given an amplitude distribution P(V < L) of a 
complex fading signal Ve’*, the mean values and the higher moments 


* The probability density function F(z, y) must be 20 for any x and y; therefore, 
the symmetric part of F(z, y) must be 2 |antisymmetric part of F(x, y)| when the 
antisymmetric part is negative. 
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of x and y are not unique. (However, the moments of the amplitude 
V are unique.) 

Therefore, specifying a Rayleigh distribution for the amplitude V 
does not necessarily imply that the mean values of x and y are zero. 
Physically if there is a direct path between the transmitter and the 
receiver of a radio link, then the mean values of x and y may not be zero. 
However, the results of this section show that the mere nonzero means 
of x and y do not necessarily exclude the Rayleigh distribution for the 
amplitude of the fading signal. 
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List of Symbols and Their Definitions 


f(a, B) 
F(a, y) 


te 
g(R) 


g.(R) 
H(a, 8) 
An, m(a, 8) 


H(—1, 0) 
H(—1, 0) 


I,(~) 
L 


m 
n 
N(L) 
0(L”) 


PV SL) 


The constant amplitude of the echo in the one-echo 
model. 

The zero order term of the Taylor Series expansion 
of ¥.(L) defined in equation (81). 

The coefficient of the nth order term of the Taylor 
series expansion V,(L) defined in equation (81). 

= byob. — bi as defined in equation (122). 

Defined by equation (121). 

The coefficient of the fourth order term L* of P(V < L) 
in equations (161) and (168). 

Defined by equation (30). 

The coefficient of the power series representation of 
P(V sS L) defined in equations (31) and (84). 

Joint probability density function of a and £. 

Joint probability density function of the real part x 
and the imaginary part y of the fading signal. 

The frequency of the sine wave. 

The probability density function of the amplitude R of 
the interfering vector. 

The nth order derivative of g(R). 

The smooth part of f(a, 8) as defined in equation (23). 
The partial derivative of H(a, 8) as defined in equa- 


tion (29). 
The value of H(a, 8) at (a = —1, 8 = 0). 
The average value of H(a, 8) at (a = —1, B = O) if 


H (qa, 8) is discontinuous at this point. 

Modified Bessel function of order n. 

An arbitrarily specified signal level in the study of the 
statistics P(V < L), N(L) and ¢(L). 

An integer. 

An integer. 

Expected number of fades per unit time below the 
specified signal level LD. 

A symbol to denote a function which goes to zero at a 
rate equal to or faster than L” as L — 0 where 7n > 0. 
Probability that the amplitude V of a fading signal 
fades below a specified signal level L. 


3258 
p(L) 
p2(V) 


pV, V) 
pilV | L) 


Vei* 
V wae 
V tad (t) 


V ver 


V.(L) 
W (@) 


W,,(8) 
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Probability density of the amplitude V at the specified 
signal level D. 

Probability density function of V. p.(V) = p(V). 
Joint probability density function of V and V. 
Conditional probability density of V under the condi- 
tion V = DL. 

The constant magnitude of the sine wave in the 
specialized model of a sine wave plus a Gaussian noise. 
In this paper, Q = 1 because of the normalization of 
the signal level. 

Joint probability density function of R and 6. 

The amplitude of the resultant interfering vector Re’’ 
of the fading signal model. 

The resultant complex interfering vector. 

An integer. 

An integer. 

A variable representing time. 

Average duration of fades below L. 

The amplitude of the envelope of a complex fading 
signal normalized to the nonfaded value V,., . 

The time derivative of the normalized amplitude V of 
the fading signal. © 

The envelope of the fading signal. 

The rms value of the time derivative of V. 

The unnormalized amplitude of the random fading 
signal. . 

The nonfaded signal level when there is no interference. 
Conditional average positive derivative of V as defined 
in equation (75). 

Probability density function of the random relative 
phase @ of one-echo model. 

The nth order derivative of the probability density 
function W (6). 

The real part of the fading signal Ve’®. 

The imaginary part of the fading signal Ve’?. 

The real part of the interfering vector Re’’. 

The imaginary part of the interfering vector Re’’. 
Gamma function. 

Defined by equation (141). 

An arbitrary constant > 0. 

The phase of the resultant interfering vector Re’’. 
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61, Defined by equation (53). 

LL (1 — p)/2 is the order of singularity of the joint prob- 
ability density function f(a, 8) at (a = —1, B = 0) as 
defined in equation (23). 

y An integer. 

pV, ¢) Joint probability density function of the amplitude V 
and the phase ¢ of the fading signal. 

a(V, ¢) A function of V and @¢ satisfying the homogeneous 
integral equation (207). 

ro) The phase of the envelope of the complex fading signal. 

Q The mean square value of an m-distributed random 
variable. 


APPENDIX B 


Amplitude Distribution In Polar Coordinates 


Since the results of P(V < L) for » = 1 cover a large class of fading 
problems and since the statistical behavior of the interfering vector is 
sometimes described by the joint probability density function q(R, 4) 
of the interfering vector, we shall also obtain the power series repre- 
sentation of P(V S JL) in terms of q(R, 6) when » = 1. By using the 
relations: 


a = R cos 8, (153) 
6 = ksin 6, (154) 


and the Jacobian relation “between f(a, 8) and q(R, 6) one can represent 
the coefficients {d.5,.} of (42) in terms of q(R, 0). This gives 


PVsbLernqi,nVtitdL+daL+-:-:-: (155) 

where 
d, = rq(l, 7) = wf(—1, 0), (156) 
dy = Zl 2) — Gol, 7) + deo(l, 7) + qoell,x)] ete, — (157) 
dan(R, 8) = =o a(R, 8). (158) 


In the deep-fade region, 
PVs Dra, wl’, (159) 
L-0 


3260 THE BELL SYSTEM TECHNICAL JOURNAL, DECEMBER 1971 


and 


p(L) oa 2rq(, 7) L. (160) 


B.1 Circular Symmetric Probability Density Function 


In this subsection, we consider a special case where the interfering 
random vector, Re’’, has a circular symmetric probability density 
function q(R, 0); i1.e., R and @ are independent and @ is uniformly dis- 
tributed in (0, 27). For this case, let g(R) be the probability density 
function of the magnitude F of the interfering vector, then 


PU SL) = $M + OL +R +--, (61 
(R, 8) = g(R) =, (162) 
C. = Flo) — wilt) + go] ete., (163) 

all) = Fe 9(R) leo (164) 


In the deep-fade region 


PW s1)—> 3g(1)L’. (165) 


These results are used in Part 3 where we discuss the relation between 
our generalized analysis and the existing theoretical work. 


APPENDIX C 


Nonanalytic H(a, B) 


In Section V of Part 1, the analysis is carried out based on the assump- 
tion that H(a, 8) can be expanded into a two-dimensional Taylor series. 
In this appendix, we shall investigate two cases where H(a, 8) cannot 
be expanded into the two-dimensional Taylor series. The objective is 
to show that from a theoretical viewpoint, the assumption of the Taylor 
series expansion of H(a, @) is not strictly necessary for the validity of 
the power laws of deep fades discussed in this paper. 


c.1 Continuous H(a, B) With Unbounded Derivatives 


If H(a, 8) is continuous at (a = —1, 8 = 0) but its first-order partial 
derivatives and/or its higher order partial derivatives are unbounded 
at (a = —1, B = 0), then H(a, B) in the neighborhood of a = —1 and 
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8 = 0 can be written as 
H(@, 8) = H(—1, 0) + O{[G + a)” + B°)"}, (166) 
where 1 > 7 > O, and 
O{{(l + a)’ + B]""} 


is a symbol to denote the component which goes to zero at a rate equal 
to or faster than that of 


(Lae) ae as. (isa) se 8] 0: 


It is obvious that H(a, 8) given by (166) cannot be expanded into 
Taylor series because the derivatives of H(a, 8) are unbounded (i.e., 
singular) at (a = —1, 6 = 0). 

Substituting equations (166) and (23) into the general formulation 
(20) for P(V s L), and carrying out the integration yields 


aH(—1, 0) 
Le 


PV<sD= Pod), (167) 


p(L) = 2nH(—1, 0)L*" + O(*"*9, (168) 


where 0(L”“*") is a symbol to denote the high-order terms which go to 
zero at a rate equal to or faster than that of L*“*" as L — 0. Since n > 0, 
then in the deep-fade region 


PV sLl)—~ 


L-0 


rH(—1, 0) lL (169) 
Be 


which is the same as equation (36). Then the discussions and conclusions 
in Sections VI, VII, and VIII of Part 1 on the power laws of deep fades 
forwy = 1,1 > yu2 1/2 andy > 1 are readily applicable to the present 
case even though the derivatives of H(a, 8) are unbounded at (a = —1, 
B = 0). 

c.2 Discontinuous H(a, 8) 


Suppose H(a, 8) and its derivatives are bounded but are discontinuous 
at 8 = 0 so that 


lim H,,,,(a, 8) # lim H,,,,,(a, 8); (170) 
bot B07 


1.e., 
FH, nla, 0°) ¥ Hy, nla, O) 
n=0,1,2,---, m=0,1,2,---. (171) 


3262 THE BELL SYSTEM TECHNICAL JOURNAL, DECEMBER 1971 


Then on each side of 8 = O, the one-sided Taylor series expansion of 
H(a, 8) is applicable. One for 8 > 0 and another for 6 < 0. Substituting 
these two Taylor series and equation (23) into the general formulation 
(20) for PV S L) and carrying out the integration, one can show that 


PV = 1) —> TALL D (172) 
10 bs 
p(L) ey 2rH(—1, 0)L**", (173) 
where 
H(—1, 0) = 3{H(—1, 0*) + A(—1, 0°)] (174) 
is the average value of the discontinuous H(a, 8) at (a2 = —1, B = 0). 


It is seen that equation (172) is also the same as equation (36) except 
for the proper interpretation of H(—1, 0) when H(a, 8) is discontinuous 
at (a2 = —1, B = 0). Therefore, the discussions and conclusions of 
Sections VI, VII, and VIII of Part 1 are also applicable to the present 
case. 


APPENDIX D 


Integration for Power Series of Amplitude Distribution 
Substituting (27) and (23) into (20) gives 


PVSL = Sa: ce: =) ae (175) 


where 





a=-14+VL2-p? 

lene = ff [1 + a)* + 6°" + a)""B" dads. (176) 
B=-L Ya=-1-VL?—$? 

From (176) it is seen that if either (n — r) or r is an odd integer, then 

I,-,,. vanishes because the integrand is antisymmetric. Therefore, 

I,-+,. does not vanish only when both (n — r) and r are even integers. 

Then let mn = 2S andr = 2». Equation (176) becomes 





@=-14+VL2—p2 —B? 

Ty5-29,20 = ~ f [1 + a)? + B7)*"(1 + a)?*~-7"8?" da dB 
B=-L Ya=-1-VL?-8? 

= 2T(S =? oe ar + 3) 2S5+2u 

> OSE OTS Do om 
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Combining (175) and (177) gives 
S=09 
PV EDS > a 
S=0 


2205 > Cop T(S — » + 3)T@ + 3) 2) = 
Es va0 (2S + 2u)T(S + 1) 25-25,2(—1, 0) | 


(178) 
S=0 

= >. dssials (179) 
S=0 


where 


__2 "SOwr(s —» + OTe + 4) 
do s+2 = (2 y! ond (28 ae 2u)T(S me 1) Hz s-2y,2.(—1, 0). (180) 


APPENDIX E 


Derivation for Expected Numbers of Fades 


Suppose 7’ = t, — t, is the time interval in which we want to find the 
expected number of times that the random fluctuating signal V(¢) 
crosses the signal level V = ZL. This interval is divided into a large 
number of smaller intervals of width At so short that each contains 
no more than one level crossing. We first consider the expected number 
of upward level crossings. The downward level crossings can be treated 
similarly. 

In an infinitesimal interval At, the conditions for an upward level 
crossing of V(t) are 


ry = Ms 0 (181) 
V(t) At > [L — Vid] > 0 (182) 


These two conditions are shown graphically in Fig. 8. On a V versus 
V plane, the region in which V and V satisfy conditions (181) and (182) 
are shown as the shaded area in Fig. 9. The integration of the joint 
probability density p(V, V) over this range will give the probability 
that V(t) will have an upward level crossing in At, 


P,.(L) ei + a _ pV, V) dV dV (183) 
~ arf" vor, | av. (184) 
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v(t) | | | 
L-v (t) 
V(t) - | 
| | | ~~-v (t) 
Lev (Ols6 “---— LEVEL L 
ore | 


| | 
— At-——>| 


t t+ At 





Fig. 8—The conditions for an upward level crossing of V(é) in an interval At. 


The expected number N,,(Z) of upward level crossings per unit time is 


Pa(L) _ 


Nup(L) _ At 


V=0 
i _ Vo, v) | dV. (185) 
v=0 V=L 
Similarly, the expected number of downward level crossings per unit 
time is 


N aown(L) =f. |V | py, V) 7 dV. (186) 


The total expected number of level crossings per unit time is 


N.(L) = N.,(L) “f= N aown(L) 


Va 
aie | V | p(¥, V) | dV. (187) 
V=-0 VeL 


Selgin les m 


Fig. 9—The region in which V and V satisfy the conditions for an upward level 
crossing. 
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The expected number of fades per unit time below V = L is 


N(L) = 3N.(L) = Nu(L) = Naown(Z) 


I 


i - Vo(¥, V) - av. (188) 


APPENDIX F 


Log Normal Behavior of M-Distributions Near the RMS Value 
Let 


7. = 90 ten (=) = 20[logie L — log. Va] (189) 
Vo 
be the signal level in dB with respect to its rms value ~/Q. Equation 
(189) implies 
Z L 
— = In (=) ; 190 
M V2 a 
L (2 ) 
—= = exp|7,], 191 
Ve Ply (191) 
aL O20 fZ\ 0b 
52 =a ie) = He as 
where 


Then the probability density functions of Z and Z are related by the 
Jacobian relation: 


na) = pt) || = 2 pe. (194) 


Substituting the m-distribution (143) and equation (191) into (194) 
yields the following probability density of Z: 


h(Z) = aes (S.)" exp [=e | 





2m” QZ 2Z | 
= MT(m) exp (m| 22 — exp (22) |\ (195) 
Substituting the following power series 
22 “1 (2Z\" 
esp (Fr) = Son Gi) ee) 
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into (195) gives 


nay = [52m ero {mf (ZF + £4 (2]}. aon 


If the signal level L is very close to its rms value +/Q, then 








z 
—-1|<Il, 198 
| V2 a 
22 L 
[22| = 2 nF] «1, (199) 
for 
|\Z|<«M/2 
= 4.3dB 
Under this condition, we have 
= 1 (2Z\" Z\ 
Sai an) <2(57) (200) 
Then h(Z) in equation (197) becomes 
MN 2m” 7m Z ‘| 
h(Z) = Ee e | exp | -2m(Z) ‘ (201) 


which is a normal distribution for Z. Substituting (190) and (201) into 
(194) yields 


p(L) = | | exp [—2m(In L — In VQ)"], (202) 


which is a log normal distribution for LZ. Therefore, within the range 
|z| << M/2 = 4.3 dB, all the m-distributions for any m 2 1/2 behave 
like a log normal distribution given by (202). 


APPENDIX G 


Nonunique Relation Between P(V S L) and f(a, 8) 


In the definition of the fading signal model discussed in Section II 
of Part 1, the four random variables V, ¢, a, and 8 are related by 


l+a= aa: 
B= Vsing 


(203) 
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Let p(V, ¢) be the joint probability density function of the amplitude 
V and the phase ¢ of the fading signal. By the Jacobian”’” of the 
transformation (203), it is easily shown that 


o(V, 4) = V0 +a)’ + fe, 8) = Vi(V cos, Vsing). (204) 


In this appendix we shall show that the relation between P(V < L) 
and p(V, ¢) is not unique. This then implies that the relation between 
P(V s L) and f(a, 8) is also nonunique because of the simple algebraic 
relation (204) between p(V, ¢) and f(a, 8). 

In the study of the amplitude distribution P(V S L) of a complex 
fading signal Ve’*, it is often assumed that V and ¢ are independent 
with @ uniformly distributed in (0, 27). However, in the study of 
interference, distortion, FM radio system, radio navigation system, etc., 
many authors”’**’**'? have investigated the distribution of the random 
phase, ¢(¢t). These results show that the distribution of phase is not 
always uniform. Furthermore, when the signal is weak (i.e., V is small), 
the phase is more likely to vary over wider range. This means the random 
variables V and ¢ are somewhat correlated. Therefore, V and ¢ generally 
can be either dependent or independent and ¢ can be either uniformly or 
nonuniformly distributed. 

By definition,”’” the probability density p(V) of V is the integration 
of p(V, ¢) over the entire range of d;i.e., 


p=2 
pV) = | a, 4) ae. (205) 


Furthermore, the cumulative amplitude distribution P(V S ZL) is the 
integration of p(V) from V = 0 to V = L. Therefore, 


pvaty=f[ [- o(V, 4) de dV. (206) 


Given a joint probability density function p(V, ¢), then ris DL) can 
be calculated by (206). 

On the other hand, given an amplitude distribution P(V S J), 
equation (206) is an integral equation to solve for p(V, ¢). An immediate 
question arising in solving the integral equation (206) is the uniqueness 
of the solution. A procedure to test the uniqueness of the solution is to 
consider the following homogeneous equation . 


V=L p=2r 
2 20 
0 i [. o(V, ¢) deaV (207) 


If the homogeneous equation (207) has a nontrivial solution, then the 
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‘solution of (206) is not unique because, given any particular solution 
pp(V, ¢) of (206), then 


p(V, 6) = pe(V, o) + colV, 9) (208) 


is also a solution of (206) where c is an arbitrary constant.* 
It is obvious that all of the following functions 


o,(V, ¢) = E(V) sin (nd) a oY) cos (n¢), n= ce 2, +3, Peat 
(209) 


are nontrivial solutions of the homogeneous equation (207) where 
£,(V) and ¢,(V) are arbitrary functions of V. Notice that the nontrivial 
solutions (209) contain both symmetric and antisymmetric functions 
of ¢. Furthermore, any arbitrary linear combination of {o,(V, $)} is 
also a solution of the homogeneous equation (207). By the experience 
of Fourier series synthesis technique, we know that the linear combina- 
tion of the set {,(V, ¢)} is able to represent a very large class of either’ 
simple or complicated functions of V and ¢. 

Therefore, given an amplitude distribution P(V S JL), the integral 
equation (206) has infinitely many different solutions. 

On the other hand, in equation (206), if one imposes the following 
two additional conditions: 


(t) V and ¢ are independent, and 
(iz) @ is uniformly distributed in (0, 27), 


then 


, 1 Tei? eas 
AV, ¢) = 5 PIV) = 5-37 PW SL) oa (210) 


is the only possible solution of (206). 
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